[Uwgarp-commits] r119 - in pkg/GARPFRM: demo vignettes

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Sat Mar 15 04:07:39 CET 2014


Author: rossbennett34
Date: 2014-03-15 04:07:15 +0100 (Sat, 15 Mar 2014)
New Revision: 119

Modified:
   pkg/GARPFRM/demo/univariate_GARCH.R
   pkg/GARPFRM/vignettes/EstimatingVolatilitiesCorrelation.Rnw
   pkg/GARPFRM/vignettes/EstimatingVolatilitiesCorrelation.pdf
Log:
Updating vignette for estimating volatility and correlation

Modified: pkg/GARPFRM/demo/univariate_GARCH.R
===================================================================
--- pkg/GARPFRM/demo/univariate_GARCH.R	2014-03-12 15:40:12 UTC (rev 118)
+++ pkg/GARPFRM/demo/univariate_GARCH.R	2014-03-15 03:07:15 UTC (rev 119)
@@ -9,7 +9,7 @@
 # Specify and fit the MSFT returns to a standard ARMA(0,0)-GARCH(1,1) model
 model0 <- uvGARCH(R, armaOrder=c(0,0))
 
-# The uvGARCH function uses the excellent rugarch packages, which has a rich
+# The uvGARCH function uses the excellent rugarch package, which has a rich
 # set of functions for analysis of fitted GARCH models. The fitted model can 
 # be extracted with the getFit function. Refer to help("uGARCHfit-class") 
 # for available all methods for the uGARCHfit object that is returned by getFit.
@@ -76,3 +76,4 @@
 model21 <- uvGARCH(R, garchOrder=c(2,1))
 getSpec(model21)
 getFit(model21)
+

Modified: pkg/GARPFRM/vignettes/EstimatingVolatilitiesCorrelation.Rnw
===================================================================
--- pkg/GARPFRM/vignettes/EstimatingVolatilitiesCorrelation.Rnw	2014-03-12 15:40:12 UTC (rev 118)
+++ pkg/GARPFRM/vignettes/EstimatingVolatilitiesCorrelation.Rnw	2014-03-15 03:07:15 UTC (rev 119)
@@ -6,7 +6,7 @@
 \begin{document}
 
 <<echo=FALSE>>=
-opts_chunk$set(tidy=FALSE, warning=FALSE, fig.width=5, fig.height=5)
+opts_chunk$set(cache=TRUE, tidy=FALSE, warning=FALSE, fig.width=5, fig.height=5)
 @
 
 
@@ -91,7 +91,7 @@
 
 Load the package and data. Unless noted otherwise, the weekly returns of Microsoft (MSFT) will be used as the asset return data.
 <<>>=
-library(GARPFRM)
+suppressPackageStartupMessages(library(GARPFRM))
 data(crsp_weekly)
 
 # Use the weekly MSFT returns
@@ -116,8 +116,10 @@
 
 Now we plot the estimated volatility from the EWMA model and the realized volatility.
 <<>>=
-plot(volEst, main="EWMA Volatility Estimate vs. Realized Volatility")
-lines(vol, col="red")
+plot(vol, main="EWMA Volatility Estimate vs. Realized Volatility")
+lines(volEst$estimate, col="red")
+legend("topright", legend=c("Realized Volatility", "EWMA Volatility Estimate"), 
+       col=c("black", "red"), lty=c(1,1), cex=0.8, bty="n")
 @
 
 The \code{estimateLambda} function estimates the value for $\lambda$ by minimizing the mean squared error between the realized volatility and the EWMA model estimated volatility.
@@ -138,17 +140,34 @@
 plot(vol, main="EWMA Volatility Estimate vs. Realized Volatility")
 # EWMA volatility estimate, lambda = 0.94
 lines(volEst$estimate, col="red")
-# EWMA volatility estimate, lambda = 0.7187
+# EWMA volatility estimate, lambda = 0.0.7359253
 lines(volEst2$estimate, col="blue")
 legend("topright", legend=c("Realized Volatility", 
                             "EWMA Volatility, lambda = 0.94",
-                            "EWMA Volatility, lambda = 0.7187"),
+                            "EWMA Volatility, lambda = 0.7359253"),
        col=c("black", "red", "blue"), lty=c(1, 1, 1), cex=0.7, bty="n")
 @
 
 Now we move to using the EWMA model to calculate the covariance between the returns of two assets. Note that we set \code{lambda=NULL} in the \code{EWMA} function. If \code{lambda = NULL}, the optimal $\lambda$ value is estimated by minimizing the mean squared error between the estimated covariance and realized covariance.
+
+The covariance between two variables, $X$ and $Y$, is defined as
+\begin{equation}
+\sigma_{XY} = E[(X - \mu_X)(Y - \mu_Y)]
+\end{equation}
+
+where
+\begin{description}
+  \item[$\mu_X$] is the sample mean of $X$
+  \item[$\mu_Y$] is the sample mean of $Y$
+\end{description}
+
+An EWMA model for the updating the covariance estimate between $X$ and $Y$ on day $n$ is
+\begin{equation}
+cov(X,Y) = \sigma_{XY,n} = \lambda \sigma_{XY,n-1} + (1 - \lambda) X_{n-1} Y_{n-1}
+\end{equation}
+
 <<>>=
-# Use the first 2 columns of the large
+# Use the first 2 columns of the large cap weekly returns
 R <- largecap_weekly[,1:2]
 initialWindow <- 52
 covEst <- EWMA(R, lambda=NULL, initialWindow, n=10, "covariance")
@@ -158,9 +177,170 @@
 
 
 In a similar fashion, we can also use the EWMA model to calculate the correlation between the returns of two assets.
+
+The correlation between two variables, $X$ and $Y$, is defined as
+\begin{equation}
+\rho_{XY} = \frac{cov(X,Y)}{\sigma_X \sigma_Y}
+\end{equation}
+
+where
+\begin{description}
+  \item[$cov(X,Y)$] is the covariance between $X$ and $Y$
+  \item[$\sigma_X$] is the standard deviation of $X$
+  \item[$\sigma_Y$] is the standard deviation of $Y$
+\end{description}
+
+The EWMA mdoel for correlation is calculated using an EWMA model for the estimated covariance between $X$ and $Y$, estimated volatility of $X$, and estimated volatility of $Y$.
 <<>>=
 corEst <- EWMA(R, lambda=NULL, initialWindow, n=10, "correlation")
 corEst
 plot(corEst, main="EWMA Estimated Correlation")
 @
+
+The previous examples demonstrated using and EWMA model to estimate the volatility of the returns of a single asset, and the correlation and volatility between the returns of two assets. Now we move to using an EWMA model to estimated the covariance and correlation of a multivariate data set.
+
+<<>>=
+# Use the first 4 columns of the largecap_weekly dataset
+R <- largecap_weekly[,1:4]
+
+# calculate the sample covariance matrix
+sample_cov <- cov(R)
+sample_cov
+
+# EWMA covariance matrix estimate
+lambda <- 0.94
+initialWindow <- 52
+covEst <- EWMA(R, lambda, initialWindow, type="covariance")
+covEst
+@
+
+TODO: GETTING UNRELIABLE CORRELATION ESTIMATES... NOT SURE WHY
+In a similar fashion, we can also use the EWMA model to estimate the correlation matrix.
+<<>>=
+# calculate the sample covariance matrix
+sample_cor <- cor(R)
+sample_cor
+
+# EWMA covariance matrix estimate
+lambda <- 0.94
+initialWindow <- 52
+corEst <- EWMA(R, lambda, initialWindow, type="correlation")
+corEst
+@
+
+
+\section{The GARCH(1,1) Model}
+We now present the generalized autoregressive conditional heteroskedasticity (GARCH) as presented by Bollerslev in 1986 as a way to estimate volatility. The general GARCH(p,q) model calculates $\sigma_n^2$ from the most recent $p$ observations of $u^2$ and the most recent $q$ estimates of $\sigma_n^2$. The GARCH(1,1) model refers to the most recent observation of $u^2$ and the most recent estimate of $\sigma_n^2$. The GARCH(1,1) is a popular model and the one we will focus on. The equation for the GARCH(1,1) model is
+
+\begin{equation}
+\sigma_n^2 = \gamma V_L + \alpha u_{n-1}^2 + \beta \sigma_{n-1}^2
+\end{equation}
+
+where:
+\begin{description}
+  \item[$\gamma$] is the weight assigned to $V_L$
+  \item[$V_L$] is the long-run average variance rate
+  \item[$alpha$] is the weight assigned to $u_{n-1}^2$
+  \item[$u_{n-1}$] is the squared returns of preiod $n-1$
+  \item[$\beta$] is the weight assigned to $\sigma_{n-1}^2$
+  \item[$\sigma_{n-1}^2$] is the estimated variance rate of period $n-1$
+\end{description}
+
+The weights must sum to 1 such that
+\begin{equation}
+\gamma + \alpha + \beta = 1
+\end{equation}
+
+It should be noted that the EWMA model discussed in the previous section is a special case of the GARCH(1,1) model where $\gamma = 0$, $\alpha = 1 - \lambda$, and $\beta = \lambda$.
+
+A more common form of the model is obtained by setting $\omega = \gamma V_L$ such that the equation for the model is
+
+\begin{equation}
+\sigma_n^2 = \omega + \alpha u_{n-1}^2 + \beta \sigma_{n-1}^2
+\end{equation}
+
+With the estimated parameters for $\omega$, $\alpha$, and $\beta$, we can calculate $\gamma$ and $V_L$ as
+\begin{eqnarray}
+\gamma &=& 1 - \alpha - \beta\\
+V_L &=& \omega / \gamma
+\end{eqnarray}
+
+A key characteristic of the GARCH(1,1) model is mean reversion, i.e. the variance rate is pulled back to the long-run average variance rate over time. In contrast, the EWMA model is not mean reverting.
+
+\subsection{Estimating GARCH(1,1) Parameters}
+Estimating the parameters for the GARCH model requires an optimization routine to maximize the likelihood. The FRM text describes an example of using a spreadsheet and a solver, e.g. the Microsoft Excel Solver. The implementation of GARCH in the GARPFRM(citation) package utilizes the rugarch(citation) package. The rugarch package uses C code for a fast and efficient algorithm for the main part of the likelihood calculation.
+
+Here we demonstrate how to specify and fit a GARCH(1,1) model using weekly returns for Microsoft.
+<<>>=
+# Use the weekly MSFT returns
+R <- largecap_weekly[,"MSFT"]
+
+# Specify and fit the MSFT returns to a standard ARMA(0,0)-GARCH(1,1) model
+# Note that the default is ARMA(1,1)-GARCH(1,1) so we only need to change
+# the ARMA order. The default arguments were chosen to be consistent with the 
+# default arguments in rugarch.
+model <- uvGARCH(R, armaOrder=c(0,0))
+
+# Get the fitted GARCH model
+fit <- getFit(model)
+
+# Get the coefficients
+coef(fit)
+
+# Show the summary results of the fit
+fit
+@
+
+\subsection{Using GARCH(1,1) to Forecast Future Volatility}
+%Recall that the GARCH(1,1) model estimates the variance rate on day $n$, $\sigma_n^2$ at the end of day $n-1$ is
+%
+%\begin{equation}
+%\sigma_n^2 = (1 - \alpha - \beta) V_L + \alpha u_{n-1}^2 + \beta \sigma_{n-1}^2
+%\end{equation}
+%
+%so that
+%\begin{equation}
+%\sigma_n^2  - V_L = \alpha ( u_{n-1}^2 - V_L) + \beta( \sigma_{n-1}^2 - V_L)
+%\end{equation}
+%
+%On day $n + t$ in the future,
+%\begin{equation}
+%\sigma_{n+t}^2  - V_L = \alpha ( u_{n+t-1}^2 - V_L) + \beta( \sigma_{n+t-1}^2 - V_L)
+%\end{equation}
+%
+%The expected value of $u_{n+t-1}^2$ is $\sigma_{n+t-1}^2$, which results in
+%\begin{equation}
+%E[ \sigma_{n+t}^2  - V_L] = (\alpha + \beta) E[ \sigma_{n+t-1}^2 - V_L ]
+%\end{equation}
+%
+%Using this equation repeatedly yields
+%\begin{equation}
+%E[ \sigma_{n+t}^2 ] = V_L + (\alpha + \beta)^t ( \sigma_{n}^2 - V_L )
+%\end{equation}
+
+<<>>=
+# n period ahead forecast
+forecast10 <- forecast(model, nAhead=10)
+forecast10
+@
+
+Plot of the forecast.
+<<>>=
+plot(forecast10, which=3)
+@
+
+
+Now we specify and fit a model with \code{outSample=100} so that we can use the last 100 data points for out of sample testing and do a rolling forecast.
+<<>>=
+model11 <- uvGARCH(R, armaOrder=c(0,0), outSample=100)
+forecast2 <- forecast(model11, nRoll=10)
+forecast2
+@
+
+Plot the rolling forecast.
+<<>>=
+plot(forecast2, which=4)
+@
+
+
 \end{document}

Modified: pkg/GARPFRM/vignettes/EstimatingVolatilitiesCorrelation.pdf
===================================================================
(Binary files differ)



More information about the Uwgarp-commits mailing list