Econometrics
Exercise 8
| Kontak | : \(\downarrow\) |
| clara.evania@student.matanauniversity.ac.id | |
| https://www.instagram.com/claraevania/ | |
| RPubs | https://rpubs.com/claradellaevania/ |
11. For this exercise use the quarterly number of arrivals to Australia from New Zealand, 1981 Q1 – 2012 Q3, from data set aus_arrivals.
library(fpp3)aus_arrivalsnewzealand = aus_arrivals %>% filter(Origin == 'NZ')a. Make a time plot of your data and describe the main features of the series.
autoplot(newzealand) + ggtitle ("Kedatangan dari Newzealand")Kedatangan Turis dari Australia ke Newzealand semakin meningkat trendnya yang terlihat dari plot diatas.
b. Create a training set that withholds the last two years of available data. Forecast the test set using an appropriate model for Holt-Winters’ multiplicative method.
aus_newzealand <- newzealand %>%
slice(1:(nrow(aus_arrivals)-(4*2)))
fit <- newzealand %>%
model(
multiplicative = ETS(Arrivals~ error("M") + trend("A") + season("M"))
)
fc <- fit %>% forecast(h = "2 years")
fc %>%
autoplot(level = NULL) +
labs(title="Data Perbandingan Perkiraan ",
y="Arrivals") +
guides(colour = guide_legend(title = "Forecast"))c. Why is multiplicative seasonality necessary here?
Multiplicative Seasonality sangat penting dikarenakan dalam kasus ini, Ukuran Pola Musiman ini berkembang dan berumbuh secara proporsional dimana dengan meningkatnya trend yang sudah digambarkan oleh plot. Selain itu dalam pemodelan dengan menggunakan Multiplicative Seasonality ini perilaku pola musiman nantinya akan diambil dan diramalkan.
d. Forecast the two-year test set using each of the following methods:
- an ETS model;
- an additive ETS model applied to a log transformed series;
- a seasonal naïve method;
- an STL decomposition applied to the log transformed data followed by an ETS model applied to the seasonally adjusted (transformed) data.
Pertimbangan <- anti_join(newzealand, aus_newzealand,
by = c("Quarter", "Origin", "Arrivals"))
fit1 <- aus_newzealand %>%
model(
'ETS model' = ETS(Arrivals),
'Additive to a Log Transformed' = ETS(log(Arrivals) ~ error("A") +
trend("A") +
season("A")),
'Seasonal Naïve Method' = SNAIVE(Arrivals),
'STL decomposition applied to the log transformed' = decomposition_model(STL(log(Arrivals)),
ETS(season_adjust))
)
fc1 <- fit1 %>%
forecast(h="2 years")
fc1 %>%
autoplot(level = NULL) +
autolayer(Pertimbangan, Arrivals) +
guides(colour=guide_legend(title="Forecast")) +
labs(title = "Pengunjung New Zealand Melatih Perbandingan Perkiraan dengan Data Sebenarnya.")e. Which method gives the best forecasts? Does it pass the residual tests?
fc1 %>% accuracy(newzealand)Berdasar pada nilai RMSE dan juga Plot, Peramalan yang paling baik adalah peramalan Model ETS dimana memiliki kesalahan yang mendekati dengan data plot yang sebenarnya
modelterbaik <- fit1 %>% select('ETS model')
modelterbaik %>%
gg_tsresiduals()augment(modelterbaik) %>%
features(.resid, ljung_box)Pada ACF Plot ataupun Uji Ijung Box , keduanya menunjukkan model ETS yang telah lulus uji residual dimana semua lag pada plot ACF berada dibawah batas dan memiliki nilai Ijung Box lebih dari 0.05
f.Compare the same four methods using time series cross-validation instead of using a training and test set. Do you come to the same conclusions?
NewZealand_CrossValidation <- newzealand %>%
slice(1:(n() - 3)) %>%
stretch_tsibble(.init = 36, .step = 3)
NewZealand_CrossValidation %>%
model(
'ETS model' = ETS(Arrivals),
'Additive to a Log Transformed' = ETS(log(Arrivals) ~ error("A") +
trend("A") +
season("A")),
'Seasonal Naïve Method' = SNAIVE(Arrivals),
'STL decomposition applied to the log transformed' = decomposition_model(STL(log(Arrivals)),
ETS(season_adjust))
) %>%
forecast(h = 3) %>%
accuracy(newzealand)12 Apply cross-validation and Compute the MSE of the resulting 44-step-ahead errors
12.1 Apply cross-validation techniques to produce 1 year ahead ETS and seasonal naïve forecasts for Portland cement production (from aus_production). Use a stretching data window with initial size of 5 years, and increment the window by one observation.
library(tsibble)
library(dplyr)
library(fpp3)
library(magrittr)
library(forecast)
library(zoo)
portland_cement_crossval <- aus_production %>%
slice(1:(n()-4)) %>%
stretch_tsibble(.init = 20, .step = 1)
port_cement_crossval_f <- portland_cement_crossval %>%
model(ETS(Cement),
SNAIVE(Cement)
) %>%
forecast(h = "1 year")12.2 Compute the MSE of the resulting 44-step-ahead errors. Comment on which forecasts are more accurate. Is this what you expected?
# port_cement_crossval_f %>%
# group_by(.id, .model) %>%
# mutate(h = row_number()) %>%
# ungroup() %>%
# accuracy(aus_production, by = c(".model", "h"))SOAL
Untuk tiga horizon pertama, hasil ETS jauh lebih baik, tetapi untuk h=4 hasilnya jauh lebih dekat. Kami berharap ETS dapat tampil lebih baik dengan rangkaian panjang seperti ini karena seharusnya tidak ada masalah dalam memperkirakan parameter dan akan menyertakan tren jika diperlukan.
13 Compare ETS(), SNAIVE() and decomposition_model(STL, ???) on the following five time series. You might need to use a Box-Cox transformation for the STL decomposition forecasts. Use a test set of three years to decide what gives the best forecasts.
Produksi bir dan batu bata dari aus_production
aus_production %>% autoplot(Beer) Tidak ada masalah varians yang terlihat, jadi kita tidak perlu transformasi box-cox.
fc <- aus_production %>%
slice(1:(nrow(aus_production)-12)) %>%
model(
'ETS' = ETS(Beer),
'Seasonal Naïve' = SNAIVE(Beer),
'STL Decomposition' = decomposition_model(STL(log(Beer)),
ETS(season_adjust))
) %>%
forecast(h = "3 years")
fc %>%
autoplot(filter_index(aus_production, "2000 Q1" ~ .),level=NULL) +
guides(colour=guide_legend(title="Perkiraan")) +
ggtitle("Produksi Beer dari aus_production")fc %>% accuracy(aus_production)Berdasarkan hasil test set, ETS dan STLM adalah yang terbaik untuk dataset ini.
tidy_bricks <- aus_production %>%
filter(!is.na(Bricks))
tidy_bricks %>% autoplot(Bricks)Tidak ada masalah varians yang terlihat, jadi kita tidak perlu transformasi box-cox.
fc <- tidy_bricks %>%
slice(1:(nrow(tidy_bricks)-12)) %>%
model(
'ETS' = ETS(Bricks),
'Seasonal Naïve' = SNAIVE(Bricks),
'STL Decomposition' = decomposition_model(STL(log(Bricks)),
ETS(season_adjust))
) %>%
forecast(h = "3 years")
fc %>%
autoplot(filter_index(tidy_bricks, "1980 Q1" ~ .),level=NULL) +
guides(colour=guide_legend(title="Forecast")) +
ggtitle("Produksi Beer dari aus_production")fc %>% accuracy(tidy_bricks)Berdasarkan hasil test set, ETS dan STLM adalah yang terbaik untuk dataset ini.
Cost of drug subsidies for diabetes (ATC2 == “A10”) and corticosteroids (ATC2 == “H02”) from PBS
subsidies <- PBS %>%
filter(ATC2 %in% c("A10", "H02")) %>%
group_by(ATC2) %>%
summarise(Cost = sum(Cost))
subsidies %>%
autoplot(Cost) +
facet_grid(vars(ATC2), scales = "free_y") +
ggtitle("Biaya subsidi obat untuk diabetes dan kortikosteroid")Ada masalah varians, jadi kita membutuhkan transformasi box-cox untuk STL. Kita perlu membagi data itu, karena lambda untuk setiap data berbeda.
# Diabetes
diabet <- subsidies %>%
filter(ATC2 %in% "A10")
lambda1 <- diabet %>%
features(Cost, features = guerrero) %>%
pull(lambda_guerrero)
fc1 <- diabet %>%
filter(Month < max(Month) - 35) %>%
model(
'ETS' = ETS(Cost),
'Seasonal Naïve' = SNAIVE(Cost),
'STL Decomposition' = decomposition_model(STL(box_cox(log(Cost), lambda1)),
ETS(season_adjust))
) %>%
forecast(h = "3 years")
# Corticosteroids
corti <- subsidies %>%
filter(ATC2 %in% "H02")
lambda2 <- corti %>%
features(Cost, features = guerrero) %>%
pull(lambda_guerrero)
fc2 <- corti %>%
filter(Month < max(Month) - 35) %>%
model(
'ETS' = ETS(Cost),
'Seasonal Naïve' = SNAIVE(Cost),
'STL Decomposition' = decomposition_model(STL(box_cox(log(Cost), lambda2)),
ETS(season_adjust))
) %>%
forecast(h = "3 years")
# Combine the data again
fc <- bind_rows(fc1,fc2)
fc %>% autoplot(subsidies, level=NULL) +
guides(colour=guide_legend(title="Forecast")) +
ggtitle("Perkiraan Selama 3 tahun")fc %>%
accuracy(subsidies) %>%
arrange(ATC2)Untuk seri A10, metode STLM tampaknya yang terbaik, sedangkan untuk seri H02, SNAIVE/STLM tampaknya yang terbaik.
Total food retailing turnover for Australia from aus_retail.
food <- aus_retail %>%
filter(Industry == "Food retailing") %>%
summarise(Turnover = sum(Turnover))
autoplot(food, Turnover) +
ggtitle('Total food retailing turnover for Australia')14 a. Use ETS() to select an appropriate model for the following series: total number of trips across Australia using tourism, the closing prices for the four stocks in gafa_stock, and the lynx series in pelt. Does it always give good forecasts?
Total number of trips across Australia using
aus_trips <- tourism %>%
summarise(Trips = sum(Trips))
fit <- aus_trips %>%
model(ETS(Trips))
report(fit)## Series: Trips
## Model: ETS(A,A,A)
## Smoothing parameters:
## alpha = 0.4495675
## beta = 0.04450178
## gamma = 0.0001000075
##
## Initial states:
## l[0] b[0] s[0] s[-1] s[-2] s[-3]
## 21689.64 -58.46946 -125.8548 -816.3416 -324.5553 1266.752
##
## sigma^2: 699901.4
##
## AIC AICc BIC
## 1436.829 1439.400 1458.267
fit %>%
forecast() %>%
autoplot(aus_trips) +
ggtitle("Forecasting of total number of trips across Australia")ETS modelâs forecasts looks reasonable.
The closing prices for the four stocks in gafa_stock
gafa_stock %>%
autoplot(Close) +
facet_grid(vars(Symbol), scales = "free_y")gafa_regular <- gafa_stock %>%
group_by(Symbol) %>%
mutate(trading_day = row_number()) %>%
ungroup() %>%
as_tsibble(index = trading_day, regular = TRUE)
fit <- gafa_regular %>%
model(
ETS(Close)
)
report(fit)fit %>%
forecast(h=50) %>%
autoplot(gafa_regular %>% group_by_key() %>% slice((n() - 100):n()))ETS modelâs forecasts looks reasonable.
pelt %>%
model(
ETS(Lynx)
) %>%
forecast(h=10) %>%
autoplot(pelt) +
ggtitle("Forecasting of PELT trading records.")14 b. Find an example where it does not work well. Can you figure out why?
As seen in the pelt data set above, ETS does not work well with cyclic data.
15Show that the point forecasts from an ETS(M,A,M) model are the same as those obtained using Holt-Wintersâ multiplicative method.
library(forecast)
?ets
library(fpp2)
library(seasonal)#ETS(M,A,M) model
usp = ets(usdeaths, model = "MAM")
#hw model
usphw = hw(usdeaths, seasonal = 'multiplicative', h=1)
accuracy(usp)## ME RMSE MAE MPE MAPE MASE
## Training set 22.42166 247.1924 193.8055 0.2077024 2.225981 0.4435586
## ACF1
## Training set 0.005044861
usdeaths.ets <- ets(usdeaths, model="MAM")
usdeaths.etsF <- forecast(usdeaths.ets, h = 1)
usdeaths.HWM <- hw(usdeaths, seasonal = 'multiplicative', h=1)
usdeaths.etsF$mean## Jan
## 1979 8233.107
usdeaths.HWM$mean## Jan
## 1979 8217.64