homework #3
Hanyue Kuang
For this homework, I use the data of monthly gas usage of Lowa. I want to compare the models of ETS, ARIMA and GARCH.
## Warning: package 'fGarch' was built under R version 3.5.2## Warning: package 'forecast' was built under R version 3.5.2## Warning: package 'data.table' was built under R version 3.5.2
ggtsdisplay(gas.ts)
# From the ACF, we can see that the data need seasonal adjusting. From the PACF, there is a significant spike at lag 1 and decaying after that.
ndiffs(gas.ts)## [1] 0nsdiffs(gas.ts)## [1] 1# For arima model, i assume the order to be (1,0,0)(2,1,1)[12]seasadj_gas <- gas.ts %>%
stl(t.window = 13, s.window = "periodic", robust = TRUE ) %>% seasadj()
autoplot(seasadj_gas)
autoplot(forecast(ets(gas.ts), h=12))
checkresiduals(ets(gas.ts))
##
## Ljung-Box test
##
## data: Residuals from ETS(M,A,A)
## Q* = 37.007, df = 5.2, p-value = 7.473e-07
##
## Model df: 16. Total lags used: 21.2summary(ets(gas.ts))## ETS(M,A,A)
##
## Call:
## ets(y = gas.ts)
##
## Smoothing parameters:
## alpha = 0.0016
## beta = 1e-04
## gamma = 0.002
##
## Initial states:
## l = 143.123
## b = -0.3216
## s = 81.8588 -7.05 -62.2069 -82.5654 -81.5803 -77.8775
## -75.3104 -41.5299 21.0762 71.55 124.6376 128.9975
##
## sigma: 0.137
##
## AIC AICc BIC
## 1062.511 1069.465 1107.789
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set -0.3302563 16.02604 12.18271 -2.508416 10.5338 0.7340863
## ACF1
## Training set 0.3976607autoplot(forecast(auto.arima(gas.ts),h=12))
checkresiduals(auto.arima(gas.ts))
##
## Ljung-Box test
##
## data: Residuals from ARIMA(1,0,0)(2,1,1)[12] with drift
## Q* = 11.379, df = 16.2, p-value = 0.7964
##
## Model df: 5. Total lags used: 21.2summary(auto.arima(gas.ts))## Series: gas.ts
## ARIMA(1,0,0)(2,1,1)[12] with drift
##
## Coefficients:
## ar1 sar1 sar2 sma1 drift
## 0.4204 -0.1900 -0.2513 -0.6711 -0.2985
## s.e. 0.0964 0.2431 0.1931 0.3289 0.0829
##
## sigma^2 estimated as 267.7: log likelihood=-400.72
## AIC=813.44 AICc=814.41 BIC=828.7
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set -0.5093539 14.99253 10.17712 -0.4178253 9.399218 0.6132369
## ACF1
## Training set 0.008890589fitgarch <- garchFit(formula = ~ garch(1, 1),data =gas.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: 106
## Recursion Init: mci
## Series Scale: 84.15247
##
## Parameter Initialization:
## Initial Parameters: $params
## Limits of Transformations: $U, $V
## Which Parameters are Fixed? $includes
## Parameter Matrix:
## U V params includes
## mu -14.81923866 14.81924 1.481924 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: 150.82488: 1.48192 0.100000 0.100000 0.800000
## 1: 150.11581: 1.49413 0.164062 0.0251148 0.813961
## 2: 150.02424: 1.49343 0.161091 0.0204345 0.810709
## 3: 149.88690: 1.48024 0.171486 1.00000e-08 0.817270
## 4: 149.87701: 1.48025 0.173346 1.00000e-08 0.819048
## 5: 149.87553: 1.47986 0.145592 1.00000e-08 0.848249
## 6: 149.87400: 1.48039 0.123437 1.00000e-08 0.871310
## 7: 149.86972: 1.48088 0.0789397 1.00000e-08 0.917257
## 8: 149.86905: 1.48095 0.0710542 1.00000e-08 0.925164
## 9: 149.86873: 1.48102 0.0634236 1.00000e-08 0.933317
## 10: 149.86870: 1.48101 0.0633557 1.00000e-08 0.933254
## 11: 149.86869: 1.48100 0.0633225 1.00000e-08 0.933341
## 12: 149.86868: 1.48076 0.0632379 1.00000e-08 0.933380
## 13: 149.86866: 1.48075 0.0630091 1.00000e-08 0.933675
## 14: 149.86864: 1.48078 0.0610899 1.00000e-08 0.935663
## 15: 149.86864: 1.48076 0.0610972 1.00000e-08 0.935670
## 16: 149.86863: 1.48075 0.0610914 1.00000e-08 0.935666
## 17: 149.86863: 1.48073 0.0610852 1.00000e-08 0.935686
## 18: 149.86863: 1.48071 0.0610439 1.00000e-08 0.935716
## 19: 149.86863: 1.48069 0.0609805 1.00000e-08 0.935795
## 20: 149.86862: 1.48052 0.0593369 1.00000e-08 0.937507
## 21: 149.86861: 1.48002 0.0598569 1.00000e-08 0.936934
## 22: 149.86861: 1.48002 0.0598710 1.00000e-08 0.936947
## 23: 149.86861: 1.48002 0.0598561 1.00000e-08 0.936959
## 24: 149.86861: 1.48007 0.0597497 1.00000e-08 0.937071
## 25: 149.86861: 1.48009 0.0597386 1.00000e-08 0.937083
## 26: 149.86861: 1.48009 0.0597368 1.00000e-08 0.937084
##
## Final Estimate of the Negative LLH:
## LLH: 619.7274 norm LLH: 5.846485
## mu omega alpha1 beta1
## 124.55313153 423.03447776 0.00000001 0.93708438
##
## R-optimhess Difference Approximated Hessian Matrix:
## mu omega alpha1 beta1
## mu -0.0151272314 0.0001297649 1.200500 8.616965e-01
## omega 0.0001297649 -0.0002311613 -1.531588 -1.561415e+00
## alpha1 1.2004996968 -1.5315884910 -10311.947645 -1.041715e+04
## beta1 0.8616965198 -1.5614154226 -10417.150815 -1.056605e+04
## attr(,"time")
## Time difference of 0.008375883 secs
##
## --- END OF TRACE ---## Warning in sqrt(diag(fit$cvar)): 产生了NaNs##
## Time to Estimate Parameters:
## Time difference of 0.03956294 secscheckresiduals(residuals(fitgarch))## Warning in modeldf.default(object): Could not find appropriate degrees of
## freedom for this model.
summary(fitgarch)##
## Title:
## GARCH Modelling
##
## Call:
## garchFit(formula = ~garch(1, 1), data = gas.ts)
##
## Mean and Variance Equation:
## data ~ garch(1, 1)
## <environment: 0x7ffca3af07d0>
## [data = gas.ts]
##
## Conditional Distribution:
## norm
##
## Coefficient(s):
## mu omega alpha1 beta1
## 1.2455e+02 4.2303e+02 1.0000e-08 9.3708e-01
##
## Std. Errors:
## based on Hessian
##
## Error Analysis:
## Estimate Std. Error t value Pr(>|t|)
## mu 1.246e+02 7.773e+00 16.02 <2e-16 ***
## omega 4.230e+02 NA NA NA
## alpha1 1.000e-08 NA NA NA
## beta1 9.371e-01 NA NA NA
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Log Likelihood:
## -619.7274 normalized: -5.846485
##
## Description:
## Sat Apr 27 00:20:21 2019 by user:
##
##
## Standardised Residuals Tests:
## Statistic p-Value
## Jarque-Bera Test R Chi^2 11.3355 0.003455624
## Shapiro-Wilk Test R W 0.8721812 4.336111e-08
## Ljung-Box Test R Q(10) 328.8864 0
## Ljung-Box Test R Q(15) 557.4682 0
## Ljung-Box Test R Q(20) 744.2641 0
## Ljung-Box Test R^2 Q(10) 87.27032 1.854072e-14
## Ljung-Box Test R^2 Q(15) 193.7214 0
## Ljung-Box Test R^2 Q(20) 203.0162 0
## LM Arch Test R TR^2 67.03657 1.143174e-09
##
## Information Criterion Statistics:
## AIC BIC SIC HQIC
## 11.76844 11.86895 11.76573 11.80918garch_fc <- predict(fitgarch, n.ahead=12)
autoplot(gas.ts, series = "Data") +
autolayer(ts(garch_fc$meanForecast, start = c(1979,11), frequency = 12),series = "Forecast") +
ggtitle("GARCH: forecast")
The GARCH model only gives the forecaste of a flat straignt line. I think I need to study further about this topic. The auto.arima has the same choice as my analysis. No matter looking from the autocorrelation and AICc, ARIMA model outshines ETS. However, GARCH model has the lowest AICcs.