Modelling Volatility

Still working on the Air Passenger data. We view our original model again.

library(forecast)
Model.Arima=Arima(AirPassengers,c(3,1,3),seasonal=c(0,1,0))
checkresiduals(Model.Arima)

## 
##  Ljung-Box test
## 
## data:  Residuals from ARIMA(3,1,3)(0,1,0)[12]
## Q* = 31.279, df = 18, p-value = 0.02672
## 
## Model df: 6.   Total lags used: 24

Test for ARCH in the residuals

We want to test for conditional heteroscedasticity in the residual. The ARCH Engle’s test is constructed based on the fact that if the residuals (defined as e[t]) are heteroscedastic, the squared residuals (e^2[t]) are autocorrelated. The first type of test is to examine whether the squares of residuals are a sequence of white noise, which is called Portmanteau Q test and similar to the Ljung-Box test on the squared residuals. The second type of test proposed by Engle (1982) is the Lagrange Multiplier test which is to fit a linear regression model for the squared residuals and examine whether the fitted model is significant. So the null hypothesis is that the squared residuals are a sequence of white noise, namely, the residuals are homoscedastic

library(aTSA)
## 
## Attaching package: 'aTSA'
## The following object is masked from 'package:forecast':
## 
##     forecast
## The following object is masked from 'package:graphics':
## 
##     identify
arch.test(arima(AirPassengers,order=c(3,1,3),seasonal=c(0,1,0)),output=TRUE)
## ARCH heteroscedasticity test for residuals 
## alternative: heteroscedastic 
## 
## Portmanteau-Q test: 
##      order    PQ p.value
## [1,]     4  6.99   0.136
## [2,]     8  8.65   0.372
## [3,]    12 11.25   0.507
## [4,]    16 14.62   0.553
## [5,]    20 22.46   0.316
## [6,]    24 23.16   0.510
## Lagrange-Multiplier test: 
##      order    LM  p.value
## [1,]     4 59.78 6.54e-13
## [2,]     8 26.67 3.81e-04
## [3,]    12 14.59 2.02e-01
## [4,]    16  7.65 9.37e-01
## [5,]    20  4.12 1.00e+00
## [6,]    24  2.96 1.00e+00

The arch test confirms the presence of arch effects. Although, the Portmanteau test in inconclusive (we fail to reject the null hypothesis because the p.value is higher than 0.05), tHe Lagrange multiplier test coefficients are significant.

Next we extract the residuals from the model Model.Arima.

Resids=Model.Arima$residuals
plot(Resids)

par(mfrow=c(1,2))
acf(Resids^2)
pacf(Resids^2)

The ACF plot shows that there might be some auto correlation on at least two lag points although the significant might be low.

Building hte Garch Model

It is quite difficult to figure out the appropriate values for p & q from the correlation plots above. So we would run the garch model fro several combinations of p ( from 1 to 9) and q (from 1 to 9) and extract the corresponding AICs from each mode.

We will select the model that gives the least AIC score and use this as our final GARCH model.

library(fGarch)
## Loading required package: timeDate
## Loading required package: timeSeries
## Loading required package: fBasics
## 
## Rmetrics Package fBasics
## Analysing Markets and calculating Basic Statistics
## Copyright (C) 2005-2014 Rmetrics Association Zurich
## Educational Software for Financial Engineering and Computational Science
## Rmetrics is free software and comes with ABSOLUTELY NO WARRANTY.
## https://www.rmetrics.org --- Mail to: info@rmetrics.org
Result=data.frame(Model="m",AIC=0)
q=0
for (i in 1:9){
  for (j in 1:9){
    q=q+1
    fit=garchFit(substitute(~garch(p,q),list(p=i,q=j)),data=Resids,trace=F)
                 
      Result=rbind(Result,data.frame(Model=paste("m-",i,"-",j),AIC=fit@fit$ics[1]))
    
  }
}
## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced

## Warning in sqrt(diag(fit$cvar)): NaNs produced
Result=Result[2:nrow(Result),]

plot(Result,col=c(1,2,3,4),pch=20)

Result[which.min(Result$AIC),]
##         Model      AIC
## AIC1 m- 1 - 2 7.428341

The best model is identified with P=2 and q=1. Next we fit the model based on these parameters and analyze it.

fit=garchFit(~garch(1,2), data=Resids,trace=F)
summary(fit)
## 
## Title:
##  GARCH Modelling 
## 
## Call:
##  garchFit(formula = ~garch(1, 2), data = Resids, trace = F) 
## 
## Mean and Variance Equation:
##  data ~ garch(1, 2)
## <environment: 0x000000001a38e200>
##  [data = Resids]
## 
## Conditional Distribution:
##  norm 
## 
## Coefficient(s):
##      mu    omega   alpha1    beta1    beta2  
## 1.59975  4.06493  0.53166  0.03029  0.51959  
## 
## Std. Errors:
##  based on Hessian 
## 
## Error Analysis:
##         Estimate  Std. Error  t value Pr(>|t|)   
## mu       1.59975     0.64437    2.483  0.01304 * 
## omega    4.06493     5.64332    0.720  0.47134   
## alpha1   0.53166     0.27855    1.909  0.05630 . 
## beta1    0.03029     0.08641    0.351  0.72594   
## beta2    0.51959     0.18577    2.797  0.00516 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Log Likelihood:
##  -529.8405    normalized:  -3.679448 
## 
## Description:
##  Tue Aug 15 09:23:20 2017 by user: tbaiy 
## 
## 
## Standardised Residuals Tests:
##                                 Statistic p-Value   
##  Jarque-Bera Test   R    Chi^2  2.305953  0.3156958 
##  Shapiro-Wilk Test  R    W      0.9805457 0.03843974
##  Ljung-Box Test     R    Q(10)  7.100603  0.7159145 
##  Ljung-Box Test     R    Q(15)  15.82434  0.3938211 
##  Ljung-Box Test     R    Q(20)  21.92675  0.3444987 
##  Ljung-Box Test     R^2  Q(10)  4.288242  0.9334085 
##  Ljung-Box Test     R^2  Q(15)  12.9225   0.6082844 
##  Ljung-Box Test     R^2  Q(20)  16.21184  0.7033983 
##  LM Arch Test       R    TR^2   8.36399   0.7560783 
## 
## Information Criterion Statistics:
##      AIC      BIC      SIC     HQIC 
## 7.428341 7.531459 7.426036 7.470243

The model confirms that the following parameters mu, alpha1 and beta2 are all significant. The Ljung-box test confirms that the residuals do not have anymore heteroscedasticity because the p-values are very large.