Forecasting Time Series

Preamble

library("e1071")     # kurtosis
library("forecast")  # Acf, Pacv
library("rugarch")   # ugarch*

Loading required package: methods

Loading required package: parallel


Attaching package: 'rugarch'

The following object is masked from 'package:stats':

    sigma

library("tseries")   # garch

source("http://ptrckprry.com/course/forecasting/code/ad.test.R") # ad.test

Simulated ARCH(2)

n <- 600
arch2.spec = ugarchspec(variance.model = list(garchOrder=c(2,0)), 
                        mean.model = list(armaOrder=c(0,0)),
                        fixed.pars=list(mu = 0, omega=0.25, alpha1=0.6, alpha2=0.35))
arch2.sim <- ugarchpath(arch2.spec, n.sim=n)
x <- drop(arch2.sim@path$seriesSim)[101:600]
time <- 1:length(x)

Time Series Plot

plot(time, x, type="l", col=2)

ACFs, PACFs

Original Series

par(mfrow=c(1,2))
Acf(x)
Pacf(x)

Squares

par(mfrow=c(1,2))
Acf(x^2)
Pacf(x^2)

Fitting an ARCH(q) with `tseries` Package

# using the garch function from library("tseries")
fit1 <- garch(x, c(0,1), trace=FALSE) # ARCH(1)
fit2 <- garch(x, c(0,2), trace=FALSE) # ARCH(2)
fit3 <- garch(x, c(0,3), trace=FALSE) # ARCH(3)
fit4 <- garch(x, c(0,4), trace=FALSE) # ARCH(4)
fit5 <- garch(x, c(0,5), trace=FALSE) # ARCH(5)

AICc with `tseries` Package

N <- length(x)
loglik0 <- -0.5 * N * (1 + log(2 * pi * mean(x^2)))
loglik1 <- logLik(fit1)
loglik2 <- logLik(fit2)
loglik3 <- logLik(fit3)
loglik4 <- logLik(fit4)
loglik5 <- logLik(fit5)

loglik <- c(loglik0, loglik1, loglik2, loglik3, loglik4, loglik5)
q <- c(0, 1, 2, 3, 4, 5)
k <- q + 1
aicc <- -2 * loglik  + 2 * k * N / (N - k - 1)

print(data.frame(q, loglik, aicc))

  q    loglik     aicc
1 0 -753.1916 1508.391
2 1 -662.8530 1329.730
3 2 -628.4104 1262.869
4 3 -665.5730 1339.227
5 4 -681.6131 1373.348
6 5 -684.6167 1381.404

Selected Model

fit.tseries <- garch(x, c(0, q[which.min(aicc)]), trace=FALSE)
summary(fit.tseries)


Call:
garch(x = x, order = c(0, q[which.min(aicc)]), trace = FALSE)

Model:
GARCH(0,2)

Residuals:
     Min       1Q   Median       3Q      Max 
-3.09809 -0.65873  0.03916  0.73043  3.16028 

Coefficient(s):
    Estimate  Std. Error  t value Pr(>|t|)    
a0   0.21804     0.03653    5.968 2.40e-09 ***
a1   0.56991     0.10664    5.344 9.08e-08 ***
a2   0.36552     0.08348    4.379 1.19e-05 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Diagnostic Tests:
    Jarque Bera Test

data:  Residuals
X-squared = 0.037617, df = 2, p-value = 0.9814


    Box-Ljung test

data:  Squared.Residuals
X-squared = 0.035092, df = 1, p-value = 0.8514

Fitting an ARCH(q) with `rugarch` Package

# ARCH(1)
spec1 <- ugarchspec(variance.model = list(garchOrder=c(1,0)), 
                   mean.model = list(armaOrder=c(0,0), include.mean=FALSE))
fit1 <- ugarchfit(spec1, x)
loglik1 <- fit1@fit$LLH

q <- c(0, 1, 2, 3, 4, 5)
loglik <- rep(NA, length(q))

for (i in 1:length(q)) {
    if (q[i] == 0) {
        N <- length(x)
        loglik[i] <- -0.5 * N * (1 + log(2 * pi * mean(x^2)))
    } else {
        spec <- ugarchspec(variance.model = list(garchOrder=c(q[i],0)), 
                           mean.model = list(armaOrder=c(0,0), include.mean=FALSE))
        fit <- ugarchfit(spec, x)
        loglik[i] <- likelihood(fit)
    }
}

k <- q + 1
aicc <- -2 * loglik  + 2 * k * N / (N - k - 1)

print(data.frame(q, loglik, aicc))

  q    loglik     aicc
1 0 -753.1916 1508.391
2 1 -663.8827 1331.790
3 2 -630.4511 1266.951
4 3 -631.1588 1270.398
5 4 -630.8577 1271.837
6 5 -631.1551 1274.481

Selected Model

spec <- ugarchspec(variance.model = list(garchOrder=c(2,0)), 
                   mean.model = list(armaOrder=c(0,0), include.mean=FALSE))
fit.rugarch <- ugarchfit(spec, x)
print(fit.rugarch)


*---------------------------------*
*          GARCH Model Fit        *
*---------------------------------*

Conditional Variance Dynamics   
-----------------------------------
GARCH Model : sGARCH(2,0)
Mean Model  : ARFIMA(0,0,0)
Distribution    : norm 

Optimal Parameters
------------------------------------
        Estimate  Std. Error  t value Pr(>|t|)
omega    0.21804    0.033429   6.5224    0e+00
alpha1   0.56991    0.101314   5.6252    0e+00
alpha2   0.36552    0.075402   4.8476    1e-06

Robust Standard Errors:
        Estimate  Std. Error  t value Pr(>|t|)
omega    0.21804    0.030190   7.2222    0e+00
alpha1   0.56991    0.089285   6.3830    0e+00
alpha2   0.36552    0.077966   4.6882    3e-06

LogLikelihood : -630.4511 

Information Criteria
------------------------------------
                   
Akaike       2.5338
Bayes        2.5591
Shibata      2.5337
Hannan-Quinn 2.5437

Weighted Ljung-Box Test on Standardized Residuals
------------------------------------
                        statistic p-value
Lag[1]                     0.1306  0.7178
Lag[2*(p+q)+(p+q)-1][2]    0.6741  0.6172
Lag[4*(p+q)+(p+q)-1][5]    2.9158  0.4227
d.o.f=0
H0 : No serial correlation

Weighted Ljung-Box Test on Standardized Squared Residuals
------------------------------------
                        statistic p-value
Lag[1]                    0.02225  0.8814
Lag[2*(p+q)+(p+q)-1][5]   0.25451  0.9881
Lag[4*(p+q)+(p+q)-1][9]   1.46554  0.9583
d.o.f=2

Weighted ARCH LM Tests
------------------------------------
            Statistic Shape Scale P-Value
ARCH Lag[3]   0.02946 0.500 2.000  0.8637
ARCH Lag[5]   0.50134 1.440 1.667  0.8832
ARCH Lag[7]   1.78756 2.315 1.543  0.7622

Nyblom stability test
------------------------------------
Joint Statistic:  0.9093
Individual Statistics:              
omega  0.63404
alpha1 0.09385
alpha2 0.23263

Asymptotic Critical Values (10% 5% 1%)
Joint Statistic:         0.846 1.01 1.35
Individual Statistic:    0.35 0.47 0.75

Sign Bias Test
------------------------------------
                   t-value    prob sig
Sign Bias          1.70346 0.08911   *
Negative Sign Bias 0.06892 0.94508    
Positive Sign Bias 0.76699 0.44345    
Joint Effect       3.76187 0.28835    


Adjusted Pearson Goodness-of-Fit Test:
------------------------------------
  group statistic p-value(g-1)
1    20     15.92      0.66261
2    30     40.12      0.08199
3    40     36.32      0.59277
4    50     57.40      0.19191


Elapsed time : 0.05638599

Conditional Variance

# Using tseries package
h.tseries <- fit.tseries$fit[,1]^2

# Using rugarch package
h.rugarch <- as.numeric(sigma(fit.rugarch)^2)

time <- 1:length(x)
par(mfrow=c(1,2))

plot(time, x, type="l", col=2)
plot(time, h.rugarch, type="l", col=3)

ARCH Residuals

# Using tseries package:
resid.tseries <- residuals(fit.tseries)

# Using rugarch package:
resid.rugarch <- as.numeric(residuals(fit.rugarch, standardize=TRUE))

# Note: resid.tseries is equal to (x / sqrt(h.tseries))
#       resid.rugarch is equal to (x / sqrt(h.rugarch))

Diagnostic Checking for ARCH

resid.arch <- resid.rugarch
plot(resid.arch, type="l", col=2)

Check White Noise

ACF, PACF

par(mfrow=c(1,2))
Acf(resid.arch)
Pacf(resid.arch)

Ljung-Box

Box.test(resid.arch, lag=1, type="Ljung-Box", fitdf=0)


    Box-Ljung test

data:  resid.arch
X-squared = 0.13064, df = 1, p-value = 0.7178

Box.test(resid.arch, lag=12, type="Ljung-Box", fitdf=0)


    Box-Ljung test

data:  resid.arch
X-squared = 14.046, df = 12, p-value = 0.2978

Box.test(resid.arch, lag=24, type="Ljung-Box", fitdf=0)


    Box-Ljung test

data:  resid.arch
X-squared = 24.608, df = 24, p-value = 0.4273

Box.test(resid.arch, lag=36, type="Ljung-Box", fitdf=0)


    Box-Ljung test

data:  resid.arch
X-squared = 37.091, df = 36, p-value = 0.4185

Check Gaussianity

Normal QQ Plot

qqnorm(resid.arch, col=2)
qqline(resid.arch, col=1)

Normality Statistics

kurtosis(resid.arch) # excess kurtosis, 0 for Gaussian

[1] 0.03379127

ad.test(resid.arch)


    Anderson-Darling normality test

data:  resid.arch
A = 0.19, p-value = 0.8989

Check Squared Residuals

Time Series Plot

plot(resid.arch^2, type="l", col=3)

ACF, PACF

par(mfrow=c(1,2))
Acf(resid.arch^2)
Pacf(resid.arch^2)

Ljung-Box

Box.test(resid.arch^2, lag=1, type="Ljung-Box", fitdf=0)


    Box-Ljung test

data:  resid.arch^2
X-squared = 0.022249, df = 1, p-value = 0.8814

Box.test(resid.arch^2, lag=12, type="Ljung-Box", fitdf=0)


    Box-Ljung test

data:  resid.arch^2
X-squared = 6.3393, df = 12, p-value = 0.898

Box.test(resid.arch^2, lag=24, type="Ljung-Box", fitdf=0)


    Box-Ljung test

data:  resid.arch^2
X-squared = 20.328, df = 24, p-value = 0.678

Box.test(resid.arch^2, lag=36, type="Ljung-Box", fitdf=0)


    Box-Ljung test

data:  resid.arch^2
X-squared = 37.893, df = 36, p-value = 0.383

Forecasting Time Series

Fitting and Diagnostic Checking for ARCH Models

Patrick Perry

April 11, 2017

Preamble

Simulated ARCH(2)

Time Series Plot

ACFs, PACFs

Fitting an ARCH(q) with `tseries` Package

AICc with `tseries` Package

Selected Model

Fitting an ARCH(q) with `rugarch` Package

Selected Model

Conditional Variance

ARCH Residuals

Diagnostic Checking for ARCH

Check White Noise

Check Gaussianity

Check Squared Residuals

Forecasting Time Series

Fitting and Diagnostic Checking for ARCH Models

Patrick Perry

April 11, 2017

Preamble

Simulated ARCH(2)

Time Series Plot

ACFs, PACFs

Fitting an ARCH(q) with tseries Package

AICc with tseries Package

Selected Model

Fitting an ARCH(q) with rugarch Package

Selected Model

Conditional Variance

ARCH Residuals

Diagnostic Checking for ARCH

Check White Noise

Check Gaussianity

Check Squared Residuals

Fitting an ARCH(q) with `tseries` Package

AICc with `tseries` Package

Fitting an ARCH(q) with `rugarch` Package