程序代写 ECON7350: Applied Econometrics for Macroeconomics and Finance

ECON7350: Applied Econometrics for Macroeconomics and Finance
Tutorial 9: Modelling Volatility – II
At the end of this tutorial you should be able to:
• construct an adequate set of models with possible TGARCH errors;

• construct an adequate set of models with possible GARCH-in-mean compo- nents;
• use R to draw inference on the presence of leverage effects from a class of TGARCH models;
• use R to draw inference on the presence of time-varying risk premia from a class of GARCH-in-mean models;
• use R to estimate and forecast volatilities based on models with TGARCH errors or GARCH-in-mean components.
Problems with Solutions
Consider the daily share prices of Merck & Co., Inc. (MRK) for the period 2 January 2001 to 23 December 2013 in the data file Merck.csv. Let {yt} denote the process of share prices. Recall that we learned how to fit ARMA-ARCH/GARCH models to data last week. We now consider extensions of these models to capture possible leverage effects and time-varying risk premia.
1. Consider the class of ARMA(pm, qm)-TGARCH(ph, qh) models for the returns rt =lnyt −lnyt−1.
(a) Construct an adequate set of models for estimating and forecasting volatilities.
Solution For this tutorial, we load the following useful packages.
library(forecast)
library(dplyr)
library(rugarch)

Next, load the data, extract the variables and plot them to get an idea on what we are working with.
mydata <- read.delim("Merck.csv", header = TRUE, sep = ",") date <- as.Date(mydata$Date, format = "%d/%m/%Y") y <- mydata$y r <- diff(log(y)) plot(date[-1], r, type = "l", xlab = "", ylab = "returns") −0.3 −0.2 −0.1 0.0 0.1 2002 2004 2006 2008 2010 2012 2014 To construct and adquate set of models, we use a similar approach to the one developed in Week 8 for the class of basic GARCH models. The additional dimension considered here is whether the model includes the “threshold term’’. Note that the rugarch package has a different naming convention for threshold GARCH models. Specifically, what we refer to as TGARCH, is called the gjr- GARCH in the package, whereas TGARCH denotes a slightly different model. This is not uncommon, so it is always important to check the documentation of any particular software / package carefully! submods <- c("sGARCH", "gjrGARCH") ARMA_TGARCH_est <- list() ic_arma_tgarch <- matrix( nrow = 4 ˆ 2 * 2 ˆ 3, ncol = 7 ) colnames(ic_arma_tgarch) <- c("pm", "qm", "ph", "qh", i <- 0; t0 <- proc.time() "thresh", "aic", "bic") for (pm in 0:3) { for (qm in 0:3) { for (ph in 0:1) { for (qh in 0:1) { for (thresh in 0:1) { i <- i + 1 ic_arma_tgarch[i, 1:5] <- c(pm, qm, ph, qh, thresh) if (ph == 0 && qh == 0) { # no such thing as a homoscedastic threshold ARMA! if (!thresh) { # for models with constant variance, the ugarchspec # and ugarchfit functions do not work well; # instead, the documentation advises to use # arfimaspec and arfimafit try(silent = T, expr = ARMA_TGARCH_mod <- arfimaspec( mean.model = list(armaOrder = c(pm, qm))) ARMA_TGARCH_est[[i]] <- arfimafit( ARMA_TGARCH_mod, r) ic_arma_tgarch[i,6:7] <- infocriteria( ARMA_TGARCH_est[[i]])[1:2] try(silent = T, expr = ARMA_TGARCH_mod <- ugarchspec( mean.model = list(armaOrder = c(pm, qm)), variance.model = list(model = submods[1 + thresh], garchOrder = c(ph, qh))) ARMA_TGARCH_est[[i]] <- ugarchfit(ARMA_TGARCH_mod, r, solver = 'hybrid') ic_arma_tgarch[i,6:7] <- infocriteria( ARMA_TGARCH_est[[i]])[1:2] cat("\n", proc.time() - t0, "\n") ## 156.89 0.21 160.81 NA NA ic_aic_arma_tgarch <- ic_arma_tgarch[ order(ic_arma_tgarch[,6]),][1:20,] ic_bic_arma_tgarch <- ic_arma_tgarch[ order(ic_arma_tgarch[,7]),][1:20,] ic_int_arma_tgarch <- intersect(as.data.frame(ic_aic_arma_tgarch), as.data.frame(ic_bic_arma_tgarch)) adq_set_arma_tgarch <- as.matrix(arrange(as.data.frame( ic_int_arma_tgarch), pm, qm, ph, qh, thresh)) adq_idx_arma_tgarch <- match(data.frame( t(adq_set_arma_tgarch[, 1:5])), data.frame(t(ic_arma_tgarch[, 1:5]))) Note that we have selected the entire intersecting set of models. The next step is to have a look at residual autocorrelation. As with the basic GARCH, one must be careful in interpreting the Ljung-Box test with conditional heteroscedasticity. Indeed, we will avoid it altogether and instead just examine the SACF. nmods <- length(adq_idx_arma_tgarch) sacf_tgarch <- matrix(nrow = nmods, ncol = 15) colnames(sacf_tgarch) <- c("pm", "qm", "ph", "qh", "thresh", 1:10) for (i in 1:nmods) sacf_tgarch[i,1:5] <- adq_set_arma_tgarch[i,1:5] sacf_tgarch[i,6:15] <- lag = 10, plot = F)$acf[2:11] All the residuals look relatively small, so we do not eliminate any specifications 4 from the existing set of models. Hence, this is our adequate set. (b) Draw inference on historic volatilities. Solution Plotting estimated volatilities is the same as for the basic GARCH (note that we still do not get confidence intervals with the rugarch package). title_tgarch <- rep(NA, times = nmods) for (i in 1:nmods) title_tgarch[i] <- paste("ARMA(", as.character(adq_set_arma_tgarch[i, 1]), ",", as.character(adq_set_arma_tgarch[i, 2]), c("", "T")[1 + adq_set_arma_tgarch[i,5]], as.character(adq_set_arma_tgarch[i, 3]), ",", as.character(adq_set_arma_tgarch[i, 4]), ")", plot(date[-1], sqrt( type = "l", xlab = "", ylab = "volatilities", ylim = c(0, 0.08), main = title_tgarch[i]) ARMA(0,0)−GARCH(1,1) 2002 2004 2006 2008 2010 2012 2014 ARMA(0,0)−TGARCH(1,1) volatilities volatilities 0.00 0.02 0.04 0.06 0.08 0.00 0.02 0.04 0.06 0.08 2002 2004 2006 2008 2010 2012 2014 ARMA(0,1)−GARCH(1,1) 2002 2004 2006 2008 2010 2012 2014 ARMA(0,1)−TGARCH(1,1) volatilities volatilities 0.00 0.02 0.04 0.06 0.08 0.00 0.02 0.04 0.06 0.08 2002 2004 2006 2008 2010 2012 2014 ARMA(0,2)−TGARCH(1,1) 2002 2004 2006 2008 2010 2012 2014 ARMA(0,3)−GARCH(1,1) volatilities volatilities 0.00 0.02 0.04 0.06 0.08 0.00 0.02 0.04 0.06 0.08 2002 2004 2006 2008 2010 2012 2014 ARMA(0,3)−TGARCH(1,1) 2002 2004 2006 2008 2010 2012 2014 ARMA(1,0)−GARCH(1,1) volatilities volatilities 0.00 0.02 0.04 0.06 0.08 0.00 0.02 0.04 0.06 0.08 2002 2004 2006 2008 2010 2012 2014 ARMA(1,0)−TGARCH(1,1) 2002 2004 2006 2008 2010 2012 2014 ARMA(1,1)−GARCH(1,1) volatilities volatilities 0.00 0.02 0.04 0.06 0.08 0.00 0.02 0.04 0.06 0.08 2002 2004 2006 2008 2010 2012 2014 ARMA(1,2)−GARCH(1,1) 2002 2004 2006 2008 2010 2012 2014 ARMA(1,2)−TGARCH(1,1) volatilities volatilities 0.00 0.02 0.04 0.06 0.08 0.00 0.02 0.04 0.06 0.08 2002 2004 2006 2008 2010 2012 2014 ARMA(1,3)−TGARCH(1,1) 2002 2004 2006 2008 2010 2012 2014 ARMA(3,0)−TGARCH(1,1) volatilities volatilities 0.00 0.02 0.04 0.06 0.08 0.00 0.02 0.04 0.06 0.08 2002 2004 2006 2008 2010 2012 2014 ARMA(3,1)−TGARCH(1,1) volatilities 0.00 0.02 0.04 0.06 0.08 2002 2004 2006 2008 2010 2012 2014 All volatility estimates generally look very similar across the specifications in the adequate set. TGARCH vs GARCH specifications generally do not produce noticeable differences for fixed pm, qm, ph and qh. However, larger ARMA lags (pm and qm) appear to be associated with larger “spikes’’ in volatilities when they do occur. (c) Draw inference on the existence of leverage effects. Solution To test for leverage effects, we focus on on the threshold parameter in models where the threshold is actually included. Note that the rugarch labels the threshold parameter γ, whereas we used the notation λ in the videos. for (i in 1:nmods) { if (adq_set_arma_tgarch[i, 5] == 1) { lambda_est <- ARMA_TGARCH_est[[ lambda_tvl <- ARMA_TGARCH_est[[ cat(paste0(title_tgarch[i], ": lambda_hat = ", round(lambda_est, 2), ", t-value = ", round(lambda_tvl, 2)), ## ARMA(0,0)-TGARCH(1,1): lambda_hat = 0.02, t-value = 4045.7 ## ARMA(0,1)-TGARCH(1,1): lambda_hat = 0.01, t-value = 15801.83 ## ARMA(0,2)-TGARCH(1,1): lambda_hat = 0.01, t-value = 22678.58 ## ARMA(0,3)-TGARCH(1,1): lambda_hat = 0.01, t-value = 50575.46 ## ARMA(1,0)-TGARCH(1,1): lambda_hat = 0.01, t-value = 16991.85 ## ARMA(1,2)-TGARCH(1,1): lambda_hat = 0.01, t-value = 3315.97 ## ARMA(1,3)-TGARCH(1,1): lambda_hat = 0.03, t-value = 4202.3 ## ARMA(3,0)-TGARCH(1,1): lambda_hat = 0.01, t-value = 4184.12 ## ARMA(3,1)-TGARCH(1,1): lambda_hat = 0.01, t-value = 11731.98 The null H0 : λ = 0 is easily rejected in favour of H1 : λ > 0 for all specifications at a very low significance level. Hence, we may infer that there is substantial evidence of “leverage effects’ ’.
It is worth pointing out that we can justify this conclusion based on the hypothesis test in TGARCH specifications only. The fact that some GARCH specifications are also included in the adequate is irrelevant for this purpose. The latter simply means that we cannot eliminate basic GARCH models on the basis of fit vs parsimony trade-off along with risk of residual autocorrelation.
(d) Forecast volatility for the four trading days past the end of the sample.
Solution Volatility forecasts are generated using the ugarchboot function, which also provides partial predictive intervals—i.e., that account for uncertainty related to unknown future shocks, but not estimation uncertainty. In fact, ugarchboot can also be instructed to compute predictive intervals that account for estimation uncertainty by specifying option method = “Full”, but this is computationally very time consuming, so we do not attempt it in this exercise.
Also, note that rugarch provides a special version of the plot function that will plot the volatility forecasts along with predictive intervals. Unfortunately, the functionality of this customised function is limited: it does not allow the user to modify the plot title, axis limits, etc.
for (i in 1:nmods) {
plot(ugarchboot(ARMA_TGARCH_est[[adq_idx_arma_tgarch[i]]],

method = “Partial”), which = 3)
Sigma Forecast
with Bootstrap Error Bands
(q: 5%, 25%, 75%, 95%)
forecast bootstrapped
2 4 6 8 10
Sigma Forecast
with Bootstrap Error Bands
(q: 5%, 25%, 75%, 95%)
forecast bootstrapped
2 4 6 8 10
0.008 0.010
0.008 0.009 0.010 0.011 0.012
GARCH model : gjrGARCH
GARCH model : sGARCH

Sigma Forecast
with Bootstrap Error Bands
(q: 5%, 25%, 75%, 95%)
forecast bootstrapped
2 4 6 8 10
Sigma Forecast
with Bootstrap Error Bands
(q: 5%, 25%, 75%, 95%)
forecast bootstrapped
2 4 6 8 10
0.008 0.010
0.008 0.010 0.012
GARCH model : gjrGARCH
GARCH model : sGARCH

Sigma Forecast
with Bootstrap Error Bands
(q: 5%, 25%, 75%, 95%)
forecast bootstrapped
2 4 6 8 10
Sigma Forecast
with Bootstrap Error Bands
(q: 5%, 25%, 75%, 95%)
forecast bootstrapped
2 4 6 8 10
0.008 0.010
0.008 0.009 0.010 0.011
GARCH model : sGARCH
GARCH model : gjrGARCH

Sigma Forecast
with Bootstrap Error Bands
(q: 5%, 25%, 75%, 95%)
forecast bootstrapped
2 4 6 8 10
Sigma Forecast
with Bootstrap Error Bands
(q: 5%, 25%, 75%, 95%)
forecast bootstrapped
2 4 6 8 10
0.008 0.009 0.010 0.011
0.008 0.010 0.012
GARCH model : sGARCH
GARCH model : gjrGARCH

Sigma Forecast
with Bootstrap Error Bands
(q: 5%, 25%, 75%, 95%)
forecast bootstrapped
2 4 6 8 10
Sigma Forecast
with Bootstrap Error Bands
(q: 5%, 25%, 75%, 95%)
forecast bootstrapped
2 4 6 8 10
0.008 0.009 0.010 0.011
0.008 0.010 0.012
GARCH model : gjrGARCH
GARCH model : sGARCH

Sigma Forecast
with Bootstrap Error Bands
(q: 5%, 25%, 75%, 95%)
forecast bootstrapped
2 4 6 8 10
Sigma Forecast
with Bootstrap Error Bands
(q: 5%, 25%, 75%, 95%)
forecast bootstrapped
2 4 6 8 10
0.008 0.009 0.010 0.011
0.007 0.009 0.011 0.013
GARCH model : gjrGARCH
GARCH model : gjrGARCH

Sigma Forecast
with Bootstrap Error Bands
(q: 5%, 25%, 75%, 95%)
forecast bootstrapped
2 4 6 8 10
There are subtle differences between the volatility forecasts generated by different GARCH and TGARCH specifications, but all agree that volatilities across the forecast horizon will remain well below the levels estimated around 2005 and then again in 2008/2009.
2. Consider a class of GARCH-M models for returns rt, with one ARCH lag and one GARCH lag.
(a) Construct an adequate set of models for estimating the risk-premium.
Solution We use the same approach as with GARCH/TGARCH models; here it is some what simpler because we fix all GARCH lags to be 1. So a specification consists of ARMA lags p, q and an indicator as whether or not the GARCH-in-mean term is included.
ARMA_GARCHM_est <- list() ic_arma_garchm <- matrix( nrow = 4 ˆ 2 * 2, ncol = 5 ) colnames(ic_arma_garchm) <- c("p", "q", "m", "aic", "bic") i <- 0; t0 <- proc.time() for (p in 0:3) for (q in 0:3) 0.008 0.009 0.010 0.011 GARCH model : gjrGARCH for (m in 0:1) { i <- i + 1 ic_arma_garchm[i, 1:3] <- c(p, q, m) try(silent = T, expr = ARMA_GARCHM_mod <- ugarchspec( mean.model = list(armaOrder = c(p, q), archm = m), variance.model = list(garchOrder = c(1, 1))) ARMA_GARCHM_est[[i]] <- ugarchfit(ARMA_GARCHM_mod, r, solver = 'hybrid') ic_arma_garchm[i,4:5] <- infocriteria( ARMA_GARCHM_est[[i]])[1:2] cat("\n", proc.time() - t0, "\n") ## 35.52 0.03 35.72 NA NA ic_aic_arma_garchm <- ic_arma_garchm[ order(ic_arma_garchm[,4]),][1:20,] ic_bic_arma_garchm <- ic_arma_garchm[ order(ic_arma_garchm[,5]),][1:20,] ic_int_arma_garchm <- intersect( as.data.frame(ic_aic_arma_garchm), as.data.frame(ic_bic_arma_garchm)) adq_set_arma_garchm <- as.matrix(arrange(as.data.frame( ic_int_arma_garchm), p, q, m)) data.frame(t(adq_set_arma_garchm[, 1:3])), adq_idx_arma_garchm <- match( data.frame(t(ic_arma_garchm[, 1:3]))) nmods <- length(adq_idx_arma_garchm) sacf_garchm <- matrix(nrow = nmods, ncol = 13) colnames(sacf_garchm) <- c("p", "q", "m", 1:10) for (i in 1:nmods) sacf_garchm[i,1:3] <- adq_set_arma_garchm[i,1:3] sacf_garchm[i,4:13] <- lag = 10, plot = F)$acf[2:11] (b) Draw inference on historic volatilities. Solution We use the same approach as with GARCH/TGARCH specifications in Question 1. title_garchm <- rep(NA, times = nmods) for (i in 1:nmods) title_garchm[i] <- paste("ARMA(", as.character(adq_set_arma_garchm[i, 1]), ",", as.character(adq_set_arma_garchm[i, 2]), c(")-GARCH(1,1)", ")-GARCHM(1,1)") plot(date[-1], sqrt( [1 + adq_set_arma_garchm[i,3]], type = "l", xlab = "", ylab = "volatilities", ylim = c(0, 0.08), main = title_garchm[i]) ARMA(0,0)−GARCH(1,1) 2002 2004 2006 2008 2010 2012 2014 ARMA(0,0)−GARCHM(1,1) volatilities volatilities 0.00 0.02 0.04 0.06 0.08 0.00 0.02 0.04 0.06 0.08 2002 2004 2006 2008 2010 2012 2014 ARMA(0,1)−GARCH(1,1) 2002 2004 2006 2008 2010 2012 2014 ARMA(0,1)−GARCHM(1,1) volatilities volatilities 0.00 0.02 0.04 0.06 0.08 0.00 0.02 0.04 0.06 0.08 2002 2004 2006 2008 2010 2012 2014 ARMA(0,3)−GARCH(1,1) 2002 2004 2006 2008 2010 2012 2014 ARMA(1,0)−GARCH(1,1) volatilities volatilities 0.00 0.02 0.04 0.06 0.08 0.00 0.02 0.04 0.06 0.08 2002 2004 2006 2008 2010 2012 2014 ARMA(1,0)−GARCHM(1,1) 2002 2004 2006 2008 2010 2012 2014 ARMA(1,1)−GARCH(1,1) volatilities volatilities 0.00 0.02 0.04 0.06 0.08 0.00 0.02 0.04 0.06 0.08 2002 2004 2006 2008 2010 2012 2014 ARMA(1,1)−GARCHM(1,1) 2002 2004 2006 2008 2010 2012 2014 ARMA(1,2)−GARCH(1,1) volatilities volatilities 0.00 0.02 0.04 0.06 0.08 0.00 0.02 0.04 0.06 0.08 2002 2004 2006 2008 2010 2012 2014 ARMA(2,0)−GARCH(1,1) 2002 2004 2006 2008 2010 2012 2014 ARMA(2,0)−GARCHM(1,1) volatilities volatilities 0.00 0.02 0.04 0.06 0.08 0.00 0.02 0.04 0.06 0.08 2002 2004 2006 2008 2010 2012 2014 ARMA(2,1)−GARCH(1,1) 2002 2004 2006 2008 2010 2012 2014 ARMA(3,0)−GARCHM(1,1) volatilities volatilities 0.00 0.02 0.04 0.06 0.08 0.00 0.02 0.04 0.06 0.08 2002 2004 2006 2008 2010 2012 2014 ARMA(3,2)−GARCH(1,1) volatilities 0.00 0.02 0.04 0.06 0.08 2002 2004 2006 2008 2010 2012 2014 All volatility estimates generally look very similar. Including the GARCH-in-mean term appears to produce only a minor difference in estimated volatilities for fixed p, q. The periods of high volatility match up to those estimated using GARCH and TGARCH models. (c) Draw inference on the existence of a time-varying risk premium. Solution To test for time-varying risk premia, we focus on on the coefficient of the GARCH-in-mean term in models where the GARCH-in-mean is actually included. This is δ in Engle, et al. (1987), but rugarch refers to it as “archm’ ’. for (i in 1:nmods) { if (adq_set_arma_garchm[i, 3] == 1) { archm_est <- ARMA_GARCHM_est[[ archm_tvl <- ARMA_GARCHM_est[[ cat(paste0(title_garchm[i], ": archm_hat = ", round(archm_est, 2), ", t-value = ", 32 round(archm_tvl, 2)), ## ARMA(0,0)-GARCHM(1,1): archm_hat = 0.05, t-value = 4.9 ## ARMA(0,1)-GARCHM(1,1): archm_hat = 0.05, t-value = 2812.95 ## ARMA(1,0)-GARCHM(1,1): archm_hat = 0.06, t-value = 2812.95 ## ARMA(1,1)-GARCHM(1,1): archm_hat = 0.05, t-value = 4.85 ## ARMA(2,0)-GARCHM(1,1): archm_hat = 0.04, t-value = 8.48 ## ARMA(3,0)-GARCHM(1,1): archm_hat = 0.05, t-value = 3528.4 The null H0 : δ = 0 is easily rejected in favour of H1 : δ > 0 at very low significance levels for all models with a GARCH-in-mean term. Hence, we conclude that there is substantial evidence of a “time-varying risk premium”.
(d) Forecast volatility for the four trading days past the end of the sample.
Solution To forecast volatilities, we again use the ugarchboot function.
for (i in 1:nmods) {
plot(ugarchboot(ARMA_GARCHM_est[[adq_idx_arma_garchm[i]]],
method = “Partial”), which = 3)

Sigma Forecast
with Bootstrap Error Bands
(q: 5%, 25%, 75%, 95%)
forecast bootstrapped
2 4 6 8 10
Sigma Forecast
with Bootstrap Error Bands
(q: 5%, 25%, 75%, 95%)

程序代写 CS代考加微信: powcoder QQ: 1823890830 Email: powcoder@163.com

Related Posts