For this week’s discussion, I picked the data set from Bank of America on yahoo finance: https://finance.yahoo.com/quote/BAC/history?p=BAC
#Reading our dataset
#Summary of our dataset
str(BOA)
## 'data.frame': 60 obs. of 7 variables:
## $ Date : chr "2015-08-01" "2015-09-01" "2015-10-01" "2015-11-01" ...
## $ Open : num 17.9 15.9 15.5 16.9 17.5 ...
## $ High : num 18.1 16.5 17.4 18.1 17.9 ...
## $ Low : num 14.6 15.2 14.6 16.9 16.5 ...
## $ Close : num 16.4 15.6 16.8 17.4 16.8 ...
## $ Adj.Close: num 15 14.2 15.4 16 15.4 ...
## $ Volume : num 2.33e+09 1.78e+09 1.85e+09 1.44e+09 1.77e+09 ...
head(BOA)
## Date Open High Low Close Adj.Close Volume
## 1 2015-08-01 17.91 18.07 14.60 16.44 15.00615 2325559300
## 2 2015-09-01 15.95 16.48 15.25 15.58 14.22115 1777912500
## 3 2015-10-01 15.52 17.44 14.63 16.78 15.36581 1848594500
## 4 2015-11-01 16.90 18.09 16.87 17.43 15.96102 1439390500
## 5 2015-12-01 17.52 17.89 16.50 16.83 15.41159 1771496200
## 6 2016-01-01 16.45 16.59 12.94 14.14 12.98476 2648604300
tail(BOA)
## Date Open High Low Close Adj.Close Volume
## 55 2020-02-01 33.00 35.45 27.70 28.50 28.12317 1072992600
## 56 2020-03-01 28.35 29.75 17.95 21.23 20.94929 2826499000
## 57 2020-04-01 19.93 25.32 19.51 24.05 23.88343 1636882500
## 58 2020-05-01 23.38 26.17 20.10 24.12 23.95295 1447546400
## 59 2020-06-01 24.28 29.01 23.02 23.75 23.58551 1801607100
## 60 2020-07-01 24.03 24.87 22.39 24.36 24.36000 1193407900
#Creating our time series
BOA.ts = ts(BOA$Adj.Close, frequency = 12, start = c(2015,8))
autoplot(BOA.ts) +xlab("Year") + ylab("Adjusting Closing Price in $")
#Splitting our time series into training and testing
train.ts = window(BOA.ts, end=c(2020,6))
test.ts = window(BOA.ts, start=c(2015,8))
#Building ARIMA Model and ETS models
M1.ARIMA = auto.arima(train.ts)
summary(M1.ARIMA)
## Series: train.ts
## ARIMA(0,1,0)
##
## sigma^2 estimated as 4.329: log likelihood=-124.79
## AIC=251.59 AICc=251.66 BIC=253.65
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 0.1456674 2.06294 1.582365 0.3593268 7.018769 0.3124246
## ACF1
## Training set -0.009918847
#Building ETS model
M1.ETS = ets(train.ts, model = "ZZZ")
summary(M1.ETS)
## ETS(M,N,N)
##
## Call:
## ets(y = train.ts, model = "ZZZ")
##
## Smoothing parameters:
## alpha = 0.9458
##
## Initial states:
## l = 14.848
##
## sigma: 0.0909
##
## AIC AICc BIC
## 326.8948 327.3312 333.1274
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 0.1569204 2.065125 1.587724 0.397414 7.046455 0.3134828 0.03991386
#Forecasting ARIMA model
fcast.M1.ARIMA = forecast(M1.ARIMA,h=12)
plot(fcast.M1.ARIMA)
lines(test.ts, col="red")
legend("topleft",lty=1,col=c("red","blue"),c("actual values","forecast"))
#Forecasting ETS Model
fcast.M1.ETS = forecast(M1.ETS,h=12)
plot(fcast.M1.ETS)
lines(test.ts, col="red")
legend("topleft",lty=1,col=c("red","blue"),c("actual values","forecast"))
#Building GARCH model – For the first time
GBOA = garch(BOA.ts)
##
## ***** ESTIMATION WITH ANALYTICAL GRADIENT *****
##
##
## I INITIAL X(I) D(I)
##
## 1 3.486087e+01 1.000e+00
## 2 5.000000e-02 1.000e+00
## 3 5.000000e-02 1.000e+00
##
## IT NF F RELDF PRELDF RELDX STPPAR D*STEP NPRELDF
## 0 1 3.660e+02
## 1 2 2.146e+02 4.14e-01 5.48e+00 1.4e-02 2.0e+03 1.0e+00 5.50e+03
## 2 4 2.142e+02 1.79e-03 1.76e-03 7.2e-04 2.0e+00 5.0e-02 1.40e-01
## 3 5 2.136e+02 2.57e-03 2.23e-03 1.1e-03 2.2e+00 1.0e-01 2.56e-03
## 4 7 2.136e+02 6.60e-05 6.61e-05 5.9e-05 1.1e+01 5.3e-03 3.47e-04
## 5 9 2.136e+02 1.10e-04 1.10e-04 1.2e-04 2.1e+00 1.1e-02 2.66e-04
## 6 10 2.136e+02 1.31e-04 1.32e-04 2.4e-04 2.6e+00 2.1e-02 1.54e-04
## 7 12 2.136e+02 4.15e-06 4.15e-06 1.5e-05 1.1e+01 1.3e-03 2.15e-05
## 8 14 2.136e+02 6.94e-06 6.94e-06 3.1e-05 2.1e+00 2.6e-03 1.76e-05
## 9 16 2.136e+02 1.17e-06 1.17e-06 6.3e-06 1.8e+01 5.2e-04 1.07e-05
## 10 18 2.136e+02 2.13e-06 2.13e-06 1.3e-05 3.0e+00 1.0e-03 9.72e-06
## 11 20 2.136e+02 3.41e-06 3.41e-06 2.6e-05 2.0e+00 2.1e-03 7.64e-06
## 12 23 2.136e+02 5.70e-08 5.70e-08 5.4e-07 1.3e+02 4.2e-05 4.32e-06
## 13 25 2.136e+02 1.13e-07 1.13e-07 1.1e-06 1.8e+01 8.4e-05 4.47e-06
## 14 27 2.136e+02 2.23e-08 2.23e-08 2.2e-07 3.3e+02 1.7e-05 4.37e-06
## 15 29 2.136e+02 4.44e-08 4.44e-08 4.4e-07 4.2e+01 3.3e-05 4.36e-06
## 16 31 2.136e+02 8.81e-08 8.81e-08 8.8e-07 2.1e+01 6.7e-05 4.31e-06
## 17 34 2.136e+02 1.75e-09 1.75e-09 1.8e-08 4.0e+03 1.3e-06 4.23e-06
## 18 36 2.136e+02 3.50e-09 3.50e-09 3.5e-08 5.0e+02 2.7e-06 4.24e-06
## 19 38 2.136e+02 6.99e-09 6.99e-09 7.0e-08 2.5e+02 5.4e-06 4.24e-06
## 20 40 2.136e+02 1.40e-09 1.40e-09 1.4e-08 5.0e+03 1.1e-06 4.23e-06
## 21 42 2.136e+02 2.80e-10 2.80e-10 2.8e-09 2.5e+04 2.1e-07 4.23e-06
## 22 44 2.136e+02 5.59e-11 5.59e-11 5.6e-10 1.3e+05 4.3e-08 4.23e-06
## 23 46 2.136e+02 1.12e-11 1.12e-11 1.1e-10 6.3e+05 8.6e-09 4.23e-06
## 24 48 2.136e+02 2.24e-11 2.24e-11 2.2e-10 7.8e+04 1.7e-08 4.23e-06
## 25 50 2.136e+02 4.47e-12 4.47e-12 4.5e-11 1.6e+06 3.4e-09 4.23e-06
## 26 52 2.136e+02 8.94e-12 8.94e-12 9.0e-11 2.0e+05 6.8e-09 4.23e-06
## 27 54 2.136e+02 1.79e-11 1.79e-11 1.8e-10 9.8e+04 1.4e-08 4.23e-06
## 28 56 2.136e+02 3.58e-12 3.58e-12 3.6e-11 2.0e+06 2.7e-09 4.23e-06
## 29 58 2.136e+02 7.15e-13 7.15e-13 7.2e-12 9.8e+06 5.5e-10 4.23e-06
## 30 61 2.136e+02 1.46e-14 1.43e-14 1.4e-13 4.9e+08 1.1e-11 4.23e-06
## 31 63 2.136e+02 2.40e-15 2.86e-15 2.9e-14 2.4e+09 2.2e-12 4.23e-06
## 32 66 2.136e+02 2.28e-14 2.29e-14 2.3e-13 7.6e+07 1.8e-11 4.23e-06
## 33 69 2.136e+02 5.32e-16 4.58e-16 4.6e-15 1.5e+10 3.5e-13 4.23e-06
## 34 71 2.136e+02 9.32e-16 9.16e-16 9.2e-15 1.9e+09 7.0e-13 4.23e-06
## 35 72 2.136e+02 -4.68e+07 1.83e-15 1.8e-14 3.8e+09 1.4e-12 4.23e-06
##
## ***** FALSE CONVERGENCE *****
##
## FUNCTION 2.135736e+02 RELDX 1.840e-14
## FUNC. EVALS 72 GRAD. EVALS 35
## PRELDF 1.832e-15 NPRELDF 4.228e-06
##
## I FINAL X(I) D(I) G(I)
##
## 1 3.486220e+01 1.000e+00 1.401e-03
## 2 9.537147e-01 1.000e+00 1.127e-01
## 3 8.032921e-13 1.000e+00 2.551e-01
summary(GBOA)
##
## Call:
## garch(x = BOA.ts)
##
## Model:
## GARCH(1,1)
##
## Residuals:
## Min 1Q Median 3Q Max
## 0.7457 0.9471 0.9956 1.0618 1.2193
##
## Coefficient(s):
## Estimate Std. Error t value Pr(>|t|)
## a0 3.486e+01 1.026e+03 0.034 0.973
## a1 9.537e-01 6.388e+00 0.149 0.881
## b1 8.033e-13 5.417e+00 0.000 1.000
##
## Diagnostic Tests:
## Jarque Bera Test
##
## data: Residuals
## X-squared = 1.8201, df = 2, p-value = 0.4025
##
##
## Box-Ljung test
##
## data: Squared.Residuals
## X-squared = 0.18546, df = 1, p-value = 0.6667
plot(GBOA)
library(fGarch) #for GARCH model
## Loading required package: timeDate
## Loading required package: timeSeries
##
## Attaching package: 'timeSeries'
## The following object is masked from 'package:zoo':
##
## time<-
## Loading required package: fBasics
library(rugarch)
## Loading required package: parallel
##
## Attaching package: 'rugarch'
## The following object is masked from 'package:stats':
##
## sigma
library(quantmod)
## Loading required package: TTR
##
## Attaching package: 'TTR'
## The following object is masked from 'package:fBasics':
##
## volatility
## Version 0.4-0 included new data defaults. See ?getSymbols.
library(rmgarch)
##
## Attaching package: 'rmgarch'
## The following objects are masked from 'package:xts':
##
## first, last
## The following objects are masked from 'package:dplyr':
##
## first, last
M1.GARCH = garchFit(formula = ~garch(1,1), data = train.ts)
##
## Series Initialization:
## ARMA Model: arma
## Formula Mean: ~ arma(0, 0)
## GARCH Model: garch
## Formula Variance: ~ garch(1, 1)
## ARMA Order: 0 0
## Max ARMA Order: 0
## GARCH Order: 1 1
## Max GARCH Order: 1
## Maximum Order: 1
## Conditional Dist: norm
## h.start: 2
## llh.start: 1
## Length of Series: 59
## Recursion Init: mci
## Series Scale: 6.27541
##
## Parameter Initialization:
## Initial Parameters: $params
## Limits of Transformations: $U, $V
## Which Parameters are Fixed? $includes
## Parameter Matrix:
## U V params includes
## mu -37.02908711 37.02909 3.702909 TRUE
## omega 0.00000100 100.00000 0.100000 TRUE
## alpha1 0.00000001 1.00000 0.100000 TRUE
## gamma1 -0.99999999 1.00000 0.100000 FALSE
## beta1 0.00000001 1.00000 0.800000 TRUE
## delta 0.00000000 2.00000 2.000000 FALSE
## skew 0.10000000 10.00000 1.000000 FALSE
## shape 1.00000000 10.00000 4.000000 FALSE
## Index List of Parameters to be Optimized:
## mu omega alpha1 beta1
## 1 2 3 5
## Persistence: 0.9
##
##
## --- START OF TRACE ---
## Selected Algorithm: nlminb
##
## R coded nlminb Solver:
##
## 0: 79.743710: 3.70291 0.100000 0.100000 0.800000
## 1: 77.303392: 3.81192 0.0800889 0.0947001 0.777029
## 2: 72.256050: 4.09168 0.0594709 0.121489 0.756411
## 3: 69.771178: 4.18718 0.0265619 0.130810 0.739854
## 4: 69.022880: 4.29806 0.0384650 0.153382 0.741550
## 5: 68.243546: 4.29522 0.0277395 0.154416 0.737282
## 6: 68.021533: 4.29738 0.0226070 0.156210 0.735014
## 7: 67.840062: 4.31641 0.0147524 0.167994 0.727283
## 8: 67.780420: 4.36378 0.0234580 0.173487 0.723136
## 9: 67.350435: 4.31263 0.0210322 0.180042 0.716275
## 10: 67.146102: 4.32431 0.0164987 0.192082 0.705725
## 11: 66.847463: 4.32239 0.0227089 0.204733 0.696350
## 12: 66.615504: 4.31998 0.0192022 0.216613 0.684821
## 13: 66.278567: 4.32056 0.0254987 0.241816 0.663091
## 14: 65.987920: 4.31958 0.0198269 0.261223 0.635921
## 15: 65.730118: 4.32874 0.0297725 0.272677 0.605741
## 16: 65.556880: 4.31724 0.0265267 0.305272 0.597714
## 17: 65.509834: 4.33676 0.0250803 0.306047 0.594648
## 18: 65.506776: 4.33826 0.0293260 0.308206 0.590523
## 19: 65.473951: 4.33668 0.0264988 0.308433 0.589205
## 20: 65.383790: 4.34157 0.0253648 0.319690 0.566664
## 21: 65.227224: 4.32624 0.0359681 0.337160 0.520657
## 22: 64.948348: 4.35401 0.0400496 0.424073 0.463584
## 23: 64.777615: 4.34120 0.0437692 0.490271 0.383110
## 24: 64.735483: 4.34122 0.0490963 0.536541 0.328511
## 25: 64.717027: 4.33182 0.0539522 0.571903 0.290392
## 26: 64.711803: 4.32491 0.0566985 0.585517 0.276344
## 27: 64.710991: 4.32295 0.0578448 0.588564 0.272757
## 28: 64.710917: 4.32290 0.0580753 0.588092 0.272847
## 29: 64.710911: 4.32303 0.0580771 0.587537 0.273197
## 30: 64.710911: 4.32306 0.0580653 0.587303 0.273317
## 31: 64.710910: 4.32304 0.0580635 0.587242 0.273314
## 32: 64.710910: 4.32303 0.0580659 0.587254 0.273288
## 33: 64.710910: 4.32302 0.0580668 0.587261 0.273281
##
## Final Estimate of the Negative LLH:
## LLH: 173.0726 norm LLH: 2.933434
## mu omega alpha1 beta1
## 27.1287506 2.2867152 0.5872611 0.2732807
##
## R-optimhess Difference Approximated Hessian Matrix:
## mu omega alpha1 beta1
## mu -3.0948203 -0.5466344 0.8779323 2.159859
## omega -0.5466344 -1.2524038 -1.8395149 -7.040958
## alpha1 0.8779323 -1.8395149 -50.0471588 -53.502921
## beta1 2.1598592 -7.0409581 -53.5029209 -94.621066
## attr(,"time")
## Time difference of 0.005975008 secs
##
## --- END OF TRACE ---
##
##
## Time to Estimate Parameters:
## Time difference of 0.04784489 secs
## Warning: Using formula(x) is deprecated when x is a character vector of length > 1.
## Consider formula(paste(x, collapse = " ")) instead.
summary(M1.GARCH)
##
## Title:
## GARCH Modelling
##
## Call:
## garchFit(formula = ~garch(1, 1), data = train.ts)
##
## Mean and Variance Equation:
## data ~ garch(1, 1)
## <environment: 0x000000002d61f650>
## [data = train.ts]
##
## Conditional Distribution:
## norm
##
## Coefficient(s):
## mu omega alpha1 beta1
## 27.12875 2.28672 0.58726 0.27328
##
## Std. Errors:
## based on Hessian
##
## Error Analysis:
## Estimate Std. Error t value Pr(>|t|)
## mu 27.1288 0.7210 37.627 <2e-16 ***
## omega 2.2867 1.7835 1.282 0.1998
## alpha1 0.5873 0.3014 1.948 0.0514 .
## beta1 0.2733 0.3044 0.898 0.3693
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Log Likelihood:
## -173.0726 normalized: -2.933434
##
## Description:
## Wed Jul 29 21:43:23 2020 by user: student
##
##
## Standardised Residuals Tests:
## Statistic p-Value
## Jarque-Bera Test R Chi^2 4.481175 0.106396
## Shapiro-Wilk Test R W 0.924492 0.001298827
## Ljung-Box Test R Q(10) 87.95505 1.365574e-14
## Ljung-Box Test R Q(15) 96.35712 6.37268e-14
## Ljung-Box Test R Q(20) 103.9441 2.463585e-13
## Ljung-Box Test R^2 Q(10) 6.05134 0.8109315
## Ljung-Box Test R^2 Q(15) 23.84875 0.06770697
## Ljung-Box Test R^2 Q(20) 26.08246 0.1631058
## LM Arch Test R TR^2 11.88778 0.4547341
##
## Information Criterion Statistics:
## AIC BIC SIC HQIC
## 6.002461 6.143311 5.994023 6.057443
#Comparing the GARCH model’s results with the other models’ results
acc.ARIMA = accuracy(fcast.M1.ARIMA, test.ts)
acc.ETS = accuracy(fcast.M1.ETS, test.ts)
#Accuracy of ARIMA model
acc.ARIMA
## ME RMSE MAE MPE MAPE MASE
## Training set 0.1456674 2.062940 1.582365 0.3593268 7.018769 0.3124246
## Test set 0.7744870 0.774487 0.774487 3.1793389 3.179339 0.1529160
## ACF1
## Training set -0.009918847
## Test set NA
#Accuracy of ETS model
acc.ETS
## ME RMSE MAE MPE MAPE MASE
## Training set 0.1569204 2.0651246 1.5877243 0.397414 7.046455 0.3134828
## Test set 0.7551898 0.7551898 0.7551898 3.100122 3.100122 0.1491059
## ACF1
## Training set 0.03991386
## Test set NA
#It appears that in terms of AIC and BIC the GARCH model does a way better job than the ETS or ARIMA model.
#Is it enough to say that the GARCH model is the best model ? Probably not because it would be nice to see how it does graphically and how it does in term of accuracy.
#With the current info, I would still rather use an ARIMA or ETS model over a GARCH model.
#Appendix - trying to forecast GARCH model
#M1.GARCH = garchFit(formula = ~garch(1,1),data=train.ts, cond.dist="QMLE")
#ug_spec = ugarchspec(mean.model = list(armaOrder = c(1,1)))
#if we want to specify a particular GARCH model, this is the code.
#ugfit = ugarchfit(spec = ug_spec, data = train.ts)
#ugfore = ugarchforecast(ugfit, n.ahead = 20)
#summary(ugfore)
#plot(ugfore)
`