[Uwgarp-commits] r118 - in pkg/GARPFRM: . R man sandbox vignettes
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Wed Mar 12 16:40:13 CET 2014
Author: rossbennett34
Date: 2014-03-12 16:40:12 +0100 (Wed, 12 Mar 2014)
New Revision: 118
Added:
pkg/GARPFRM/R/rollFUN.R
pkg/GARPFRM/sandbox/loop_rollapply_benchmark.R
Modified:
pkg/GARPFRM/DESCRIPTION
pkg/GARPFRM/NAMESPACE
pkg/GARPFRM/R/EWMA.R
pkg/GARPFRM/man/EWMA.Rd
pkg/GARPFRM/man/estimateLambdaCor.Rd
pkg/GARPFRM/man/estimateLambdaCov.Rd
pkg/GARPFRM/man/estimateLambdaVol.Rd
pkg/GARPFRM/man/getCor.Rd
pkg/GARPFRM/man/getCov.Rd
pkg/GARPFRM/man/plot.EWMA.Rd
pkg/GARPFRM/man/realizedCor.Rd
pkg/GARPFRM/man/realizedCov.Rd
pkg/GARPFRM/man/realizedVol.Rd
pkg/GARPFRM/sandbox/ross_EWMA.R
pkg/GARPFRM/vignettes/EstimatingVolatilitiesCorrelation.Rnw
pkg/GARPFRM/vignettes/EstimatingVolatilitiesCorrelation.pdf
Log:
Updating EWMA code and docs
Modified: pkg/GARPFRM/DESCRIPTION
===================================================================
--- pkg/GARPFRM/DESCRIPTION 2014-03-11 23:04:15 UTC (rev 117)
+++ pkg/GARPFRM/DESCRIPTION 2014-03-12 15:40:12 UTC (rev 118)
@@ -21,3 +21,4 @@
'garch11.R'
'monte_carlo.R'
'efficient_frontier.R'
+ 'rollFUN.R'
Modified: pkg/GARPFRM/NAMESPACE
===================================================================
--- pkg/GARPFRM/NAMESPACE 2014-03-11 23:04:15 UTC (rev 117)
+++ pkg/GARPFRM/NAMESPACE 2014-03-12 15:40:12 UTC (rev 118)
@@ -19,7 +19,6 @@
export(getFit)
export(getSpec)
export(getStatistics)
-export(getVar)
export(hypTest)
export(minVarPortfolio)
export(monteCarlo)
@@ -32,9 +31,6 @@
export(realizedCov)
export(realizedVol)
export(tangentPortfolio)
-export(uvEWMAcor)
-export(uvEWMAcov)
-export(uvEWMAvol)
export(uvGARCH)
S3method(countViolations,xts)
S3method(fcstGarch11,DCCfit)
Modified: pkg/GARPFRM/R/EWMA.R
===================================================================
--- pkg/GARPFRM/R/EWMA.R 2014-03-11 23:04:15 UTC (rev 117)
+++ pkg/GARPFRM/R/EWMA.R 2014-03-12 15:40:12 UTC (rev 118)
@@ -16,6 +16,43 @@
#' @param n number of periods used to calculate realized volatility, covariance, or correlation.
#' @param type estimate volatility, covariance, or correlation.
#'
+#' @examples
+#' # data and parameters for EWMA estimate
+#' data(crsp_weekly)
+#' R <- largecap_weekly[, 1:2]
+#' mvR <- largecap_weekly[,1:4]
+#' lambda <- 0.94
+#' initialWindow <- 150
+#'
+#' # volatility estimate of univariate data
+#' lambda <- estimateLambdaVol(R[,1], initialWindow, n=10)
+#' vol1 <- EWMA(R[,1], lambda=NULL, initialWindow, n=10, "volatility")
+#' vol1a <- EWMA(R[,1], lambda, initialWindow, n=10, "volatility")
+#' all.equal(vol1$estimate, vol1a$estimate)
+#' vol1
+#'
+#' # covariance estimate of bivariate data
+#' lambda <- estimateLambdaCov(R, initialWindow, n=10)
+#' cov1 <- EWMA(R, lambda=NULL, initialWindow, n=10, "covariance")
+#' cov1a <- EWMA(R, lambda, initialWindow, n=10, "covariance")
+#' all.equal(cov1$estimate, cov1a$estimate)
+#' cov1
+#'
+#' # correlation estimate of bivariate data
+#' lambda <- estimateLambdaCor(R, initialWindow, n=10)
+#' cor1 <- EWMA(R, lambda=NULL, initialWindow, n=10, "correlation")
+#' cor1a <- EWMA(R, lambda, initialWindow, n=10, "correlation")
+#' all.equal(cor1$estimate, cor1a$estimate)
+#' cor1
+#'
+#' # Multivariate EWMA estimate of covariance
+#' lambda <- 0.94
+#' cov_mv <- EWMA(mvR, lambda, initialWindow, type="covariance")
+#' cov_mv
+#'
+#' # Multivariate EWMA estimate of correlation
+#' cor_mv <- EWMA(mvR, lambda, initialWindow, type="correlation")
+#' cor_mv
#' @export
EWMA <- function(R, lambda=0.94, initialWindow=10, n=10, type=c("volatility", "covariance", "correlation")){
type <- match.arg(type)
@@ -57,7 +94,7 @@
lambda <- estimateLambdaCov(R=R, initialWindow=initialWindow, n=n)
}
est <- uvEWMAcov(R=R, lambda=lambda, initialWindow=initialWindow)
- #realized <- realizedCov(R=R, n=n)
+ realized <- realizedCov(R=R, n=n)
class <- "uvEWMAcov"
}
@@ -67,7 +104,7 @@
lambda <- estimateLambdaCor(R=R, initialWindow=initialWindow, n=n)
}
est <- uvEWMAcor(R=R, lambda=lambda, initialWindow=initialWindow)
- #realized <- realizedCor(R=R, n=n)
+ realized <- realizedCor(R=R, n=n)
class <- "uvEWMAcor"
}
@@ -108,16 +145,7 @@
}
-#' EWMA Volatility Estimate
-#'
-#' EWMA model to estimate volatility
-#'
-#' @param R xts object of asset returns
-#' @param lambda smoothing parameter, must be greater than 0 or less than 1
-#' @param initialWindow initial window of observations used in estimating the
-#' initial conditions
-#'
-#' @export
+# univariate EWMA volatility estimate
uvEWMAvol <- function(R, lambda=0.94, initialWindow=10){
# function to compute EWMA volatility estimates of univariate returns
@@ -155,16 +183,7 @@
return(out)
}
-#' EWMA Covariance Estimate
-#'
-#' EWMA model to estimate covariance
-#'
-#' @param R xts object asset returns
-#' @param lambda smoothing parameter, must be greater than 0 or less than 1
-#' @param initialWindow initial window of observations used in estimating the
-#' initial conditions
-#'
-#' @export
+# univariate EWMA covariance estimate
uvEWMAcov <- function(R, lambda=0.94, initialWindow=10){
# function to compute EWMA covariance estimates
# Check for lambda between 0 and 1 & initialWindow must be greater than ncol(R)
@@ -203,16 +222,7 @@
return(out)
}
-#' EWMA Correlation Estimate
-#'
-#' EWMA model to estimate correlation
-#'
-#' @param R xts object of asset returns
-#' @param lambda smoothing parameter, must be greater than 0 or less than 1
-#' @param initialWindow initial window of observations used in estimating the
-#' initial conditions
-#'
-#' @export
+# univariate EWMA correlation estimate
uvEWMAcor <- function(R, lambda=0.94, initialWindow=10){
# function to compute EWMA correlation estimates of a bivariate dataset
# Check for lambda between 0 and 1 & initialWindow must be greater than ncol(R)
@@ -264,6 +274,7 @@
return(out)
}
+# multivariate EWMA covariance estimate
mvEWMAcov <- function(R, lambda, initialWindow){
# Separate data into a initializing window and a testing window
initialR = R[1:initialWindow,]
@@ -285,6 +296,7 @@
return(est)
}
+# multivariate EWMA correlation estimate
mvEWMAcor <- function(R, lambda, initialWindow){
cov_est <- mvEWMAcov(R=R, lambda=lambda, initialWindow=initialWindow)
est <- lapply(cov_est, cov2cor)
@@ -297,76 +309,95 @@
#'
#' Calculate realized volatility
#'
+#' Realized volatility is defined as the standard deviation of using the
+#' previous \code{n} periods.
+#'
#' @param R xts object of asset returns
#' @param n number of periods used to calculate realized volatility
#'
+#' @examples
+#' data(crsp_weekly)
+#' R <- largecap_weekly[, 1]
+#' # Calculate realized volatility
+#' realizedVolatility <- realizedVol(R[,1], 10)
+#' head(realizedVolatility)
#' @export
realizedVol <- function(R, n){
- lag(rollapply(R[,1], width=n, FUN=sd))
+ n <- as.integer(n)[1]
+ lag(rollSD(R=R[,1], width=n))
}
-cor.fun <- function(x){
- cor(x)[1,2]
-}
-
-cov.fun <- function(x){
- cov(x)[1,2]
-}
-
#' Realized Covariance
#'
#' Calculate realized covariance
#'
+#' Realized covariance is calculated using the previous \code{n} periods.
+#'
#' @param R xts object of asset returns
#' @param n number of periods used to calculate realized volatility
#'
+#' @examples
+#' data(crsp_weekly)
+#' R <- largecap_weekly[, 1:2]
+#' # Calculate realized covariance
+#' realizedCovariance <- realizedCov(R, 10)
+#' head(realizedCovariance)
#' @export
realizedCov <- function(R, n){
- lag(rollapply(R, width=n, FUN=cov.fun))
+ n <- as.integer(n)[1]
+ lag(rollCov(R=R, width=n))
}
#' Realized Correlation
#'
#' Calculate realized correlation
#'
+#' Realized correlation is calculated using the previous \code{n} periods.
+#'
#' @param R xts object of asset returns
#' @param n number of periods used to calculate realized volatility
-#'
+#' @examples
+#' data(crsp_weekly)
+#' R <- largecap_weekly[, 1:2]
+#' # Calculate realized correlation
+#' realizedCorrelation <- realizedCor(R, 10)
+#' head(realizedCorrelation)
#' @export
realizedCor <- function(R, n){
- lag(rollapply(R, width=n, FUN=cor.fun))
+ n <- as.integer(n)[1]
+ lag(rollCor(R=R, width=n))
}
# objective function for calculating lambda
-# objective is the sum of squared differences between estimated volatility and
+# objective is the mean squared error between estimated volatility and
# realized volatility
objLambdaVol <- function(lambda, R, initialWindow, n){
realized <- realizedVol(R, n)
est <- uvEWMAvol(R, lambda, initialWindow)
tmpDat <- na.omit(cbind(est, realized))
- obj <- sum((tmpDat[,1] - tmpDat[,2])^2)
+ obj <- mean((tmpDat[,1] - tmpDat[,2])^2)
return(obj)
}
# objective function for calculating lambda
-# objective is the sum of squared differences between estimated covariance and
+# objective is the mean squared error between estimated covariance and
# realized covariance
objLambdaCov <- function(lambda, R, initialWindow, n){
realized <- realizedCov(R, n)
est <- uvEWMAcov(R, lambda, initialWindow)
tmpDat <- na.omit(cbind(est, realized))
- obj <- sum((tmpDat[,1] - tmpDat[,2])^2)
+ obj <- mean((tmpDat[,1] - tmpDat[,2])^2)
return(obj)
}
# objective function for calculating lambda
-# objective is the sum of squared differences between estimated correlation and
+# objective is the mean squared error between estimated correlation and
# realized correlation
objLambdaCor <- function(lambda, R, initialWindow, n){
realized <- realizedCor(R, n)
est <- uvEWMAcor(R, lambda, initialWindow)
tmpDat <- na.omit(cbind(est, realized))
- obj <- sum((tmpDat[,1] - tmpDat[,2])^2)
+ obj <- mean((tmpDat[,1] - tmpDat[,2])^2)
return(obj)
}
@@ -374,11 +405,19 @@
#'
#' Estimate lambda for EWMA volatility estimate
#'
+#' The optimal value for lambda is calcualted by minimizing the mean squared
+#' error between the estimated volatility and realized volatility.
+#'
#' @param R xts object of asset returns
#' @param initialWindow initial window of observations used in estimating the
#' initial
#' @param n number of periods used to calculate realized volatility
#'
+#' @examples
+#' data(crsp_weekly)
+#' R <- largecap_weekly[, 1]
+#' initialWindow <- 150
+#' lambda <- estimateLambdaVol(R, initialWindow, n=10)
#' @export
estimateLambdaVol <- function(R, initialWindow=10, n=10){
opt <- optimize(objLambdaVol, interval=c(0,1), R=R,
@@ -392,11 +431,19 @@
#'
#' Estimate lambda for EWMA covariance estimate
#'
+#' The optimal value for lambda is calcualted by minimizing the mean squared
+#' error between the estimated covariance and realized covariance.
+#'
#' @param R xts object of asset returns
#' @param initialWindow initial window of observations used in estimating the
#' initial
#' @param n number of periods used to calculate realized covariance
#'
+#' @examples
+#' data(crsp_weekly)
+#' R <- largecap_weekly[, 1:2]
+#' initialWindow <- 150
+#' lambda <- estimateLambdaCov(R, initialWindow, n=10)
#' @export
estimateLambdaCov <- function(R, initialWindow=10, n=10){
opt <- optimize(objLambdaCov, interval=c(0,1), R=R,
@@ -410,11 +457,18 @@
#'
#' Estimate lambda for EWMA correlation estimate
#'
+#' The optimal value for lambda is calcualted by minimizing the mean squared
+#' error between the estimated correlation and realized correlation.
+#'
#' @param R xts object of asset returns
#' @param initialWindow initial window of observations used in estimating the
#' initial
#' @param n number of periods used to calculate realized correlation
-#'
+#' @examples
+#' data(crsp_weekly)
+#' R <- largecap_weekly[, 1:2]
+#' initialWindow <- 150
+#' lambda <- estimateLambdaCor(R, initialWindow, n=10)
#' @export
estimateLambdaCor <- function(R, initialWindow=10, n=10){
opt <- optimize(objLambdaCor, interval=c(0,1), R=R,
@@ -439,13 +493,22 @@
}
# extract the covariance between two assets from an mvEWMAcov object
+
#' EWMA Covariance
#'
-#' Extract the covariance of two assets from an \code{EWMA} object
+#' Extract the covariance of two assets from an \code{mvEWMAcov} object
#'
#' @param object an EWMA object created by \code{EWMA}
#' @param assets character vector or numeric vector. The assets can be
#' specified by name or index.
+#' @examples
+#' data(crsp_weekly)
+#' mvR <- largecap_weekly[,1:4]
+#' lambda <- 0.94
+#' initialWindow <- 150
+#' cov_mv <- EWMA(mvR, lambda, initialWindow, type="covariance")
+#' # Extract the estimated covariance between ORCL and MSFT
+#' tail(getCov(cov_mv, assets=c("ORCL", "MSFT")))
#' @export
getCov <- function(EWMA, assets){
UseMethod("getCov")
@@ -479,46 +542,54 @@
}
-# extract the variance from an mvEWMAcov object
-#' EWMA Variance
-#'
-#' Extract the Variance of an asset from an \code{EWMA} object
-#'
-#' @param object an EWMA object created by \code{EWMA}
-#' @param assets character vector or numeric vector. The assets can be
-#' specified by name or index.
-#' @export
-getVar <- function(EWMA, assets){
- UseMethod("getVar")
-}
+# # extract the variance from an mvEWMAcov object
+# # EWMA Variance
+# #
+# # Extract the Variance of an asset from an \code{mvEWMAcov} object
+# #
+# # @param object an EWMA object created by \code{EWMA}
+# # @param assets character vector or numeric vector. The assets can be
+# # specified by name or index.
+# # @export
+# getVar <- function(EWMA, assets){
+# UseMethod("getVar")
+# }
+#
+# # @method getCov mvEWMAcov
+# # @S3method getCov mvEWMAcov
+# getVar.mvEWMAcov <- function(EWMA, assets=1){
+# if(!inherits(EWMA, "mvEWMAcov")) stop("EWMA must be of class mvEWMAcov")
+#
+# cnames <- colnames(EWMA$estimate[[1]])
+#
+# # Check if asset is a character
+# if(is.character(assets[1])){
+# idx1 = grep(assets[1], cnames)
+# if(length(idx1) == 0) stop("name for asset not in EWMA object")
+# } else {
+# idx1 = assets[1]
+# }
+# out = xts(unlist(lapply(EWMA$estimate, function(x) x[idx1])), as.Date(names(EWMA$estimate)))
+# colnames(out) = cnames[idx1]
+# return(out)
+# }
-#' @method getCov mvEWMAcov
-#' @S3method getCov mvEWMAcov
-getVar.mvEWMAcov <- function(EWMA, assets=1){
- if(!inherits(EWMA, "mvEWMAcov")) stop("EWMA must be of class mvEWMAcov")
-
- cnames <- colnames(EWMA$estimate[[1]])
-
- # Check if asset is a character
- if(is.character(assets[1])){
- idx1 = grep(assets[1], cnames)
- if(length(idx1) == 0) stop("name for asset not in EWMA object")
- } else {
- idx1 = assets[1]
- }
- out = xts(unlist(lapply(EWMA$estimate, function(x) x[idx1])), as.Date(names(EWMA$estimate)))
- colnames(out) = cnames[idx1]
- return(out)
-}
-
# extract the correlation between two assets from an mvEWMAcor object
#' EWMA Correlation
#'
-#' Extract the correlation of two assets from an \code{EWMA} object
+#' Extract the correlation of two assets from an \code{mvEWMAcor} object
#'
#' @param object an EWMA object created by \code{EWMA}
#' @param assets character vector or numeric vector. The assets can be
#' specified by name or index.
+#' @examples
+#' data(crsp_weekly)
+#' mvR <- largecap_weekly[,1:4]
+#' lambda <- 0.94
+#' initialWindow <- 150
+#' cor_mv <- EWMA(mvR, lambda, initialWindow, type="correlation")
+#' # Extract the estimated correlation between ORCL and MSFT
+#' tail(getCor(cor_mv, assets=c("ORCL", "MSFT")))
#' @export
getCor <- function(EWMA, assets){
UseMethod("getCor")
@@ -551,44 +622,70 @@
return(out)
}
-# plot methods for
-# uvEWMAvol
-# uvEWMAcov
-# uvEWMAcor
-# mvEWMAcov
-# mvEWMAcor
-
#' Plot EWMA Model Estimates
#'
-#' Plot method for EWMA objects
+#' Plot method for EWMA objects.
#'
#' @param x an EWMA object
#' @param y NULL
#' @param \dots additional arguments passed to \code{plot.xts}
-#' @param assets for multivariate EWMA objects, character vector or numeric
-#' vector of assets to extract from the covariance or correlation matrix.
-#' The assets can be specified by name or index.
+#' @param assets character vector or numeric vector of assets to extract from
+#' the covariance or correlation matrix. The assets can be specified by name or
+#' index. This argument is only usd for multivariate EWMA estimates of
+#' a covariance or correlation matrix.
#' @param legendLoc location of legend. If NULL, the legend will be omitted
#' from the plot
#' @param main main title for the plot
+#' @examples
+#' # data and parameters for EWMA estimate
+#' data(crsp_weekly)
+#' R <- largecap_weekly[, 1:2]
+#' mvR <- largecap_weekly[,1:4]
+#' lambda <- 0.94
+#' initialWindow <- 150
+#'
+#' # volatility estimate of univariate data
+#' vol1 <- EWMA(R[,1], lambda, initialWindow, type="volatility")
+#' plot(vol1)
+#'
+#' # covariance estimate of bivariate data
+#' cov1 <- EWMA(R, lambda, initialWindow, type="covariance")
+#' plot(cov1)
+#'
+#' # correlation estimate of bivariate data
+#' cor1 <- EWMA(R, lambda, initialWindow, type="correlation")
+#' plot(cor1)
+#'
+#' # Multivariate EWMA estimate of covariance
+#' cov_mv <- EWMA(mvR, lambda, initialWindow, type="covariance")
+#' # These two are equivalent
+#' plot(cov_mv, assets=c("ORCL", "MSFT"))
+#' plot(cov_mv, assets=c(1, 2))
+#'
+#' # Multivariate EWMA estimate of correlation
+#' cor_mv <- EWMA(mvR, lambda, initialWindow, type="correlation")
+#' # These two are equivalent
+#' plot(cor_mv, assets=c("ORCL", "EMC"))
+#' plot(cor_mv, assets=c(1, 4))
#' @method plot EWMA
#' @S3method plot EWMA
plot.EWMA <- function(x, y=NULL, ..., assets=c(1,2), legendLoc=NULL, main="EWMA Estimate", cexLegend=0.8){
-
+ type <- x$model$type
if(inherits(x, "uvEWMAvol") | inherits(x, "uvEWMAcov") | inherits(x, "uvEWMAcor")){
- # uvEWMA has same format
+ # all uvEWMA* objects have same format
estValues <- x$estimate
} else if(inherits(x, "mvEWMAcov") | inherits(x, "mvEWMAcor")){
- # mvEWMA stuff
+ # the covariance of two assets must be extracted
if(inherits(x, "mvEWMAcov")){
estValues <- getCov(x, assets)
}
-
+ # the correlation of two assets must be extracted
if(inherits(x, "mvEWMAcor")){
estValues <- getCor(x, assets)
}
}
- plot.xts(x=estValues, ...=..., type="l", main=main)
+ # plot the values
+ plot.xts(x=estValues, ...=..., type="l", ylab=type, main=main)
if(!is.null(legendLoc)){
legend(legendLoc, legend=c("EWMA Estimate"),
lty=1, col="black", bty="n", cex=cexLegend)
Added: pkg/GARPFRM/R/rollFUN.R
===================================================================
--- pkg/GARPFRM/R/rollFUN.R (rev 0)
+++ pkg/GARPFRM/R/rollFUN.R 2014-03-12 15:40:12 UTC (rev 118)
@@ -0,0 +1,41 @@
+rollCov <- function(R, width){
+ n <- nrow(R)
+ out <- xts(vector("numeric", n), index(R))
+ for(i in width:n){
+ tmpR <- R[(i-width+1):i,]
+ out[i] <- cov(tmpR[,1], tmpR[,2])
+ }
+ # pad with leading NA
+ for(i in 1:(width-1)){
+ out[i] <- NA
+ }
+ out
+}
+
+rollCor <- function(R, width){
+ n <- nrow(R)
+ out <- xts(vector("numeric", n), index(R))
+ for(i in width:n){
+ tmpR <- R[(i-width+1):i,]
+ out[i] <- cor(tmpR[,1], tmpR[,2])
+ }
+ # pad with leading NA
+ for(i in 1:(width-1)){
+ out[i] <- NA
+ }
+ out
+}
+
+rollSD <- function(R, width){
+ n <- nrow(R)
+ out <- xts(vector("numeric", n), index(R))
+ for(i in width:n){
+ tmpR <- R[(i-width+1):i,]
+ out[i] <- sd(tmpR[,1])
+ }
+ # pad with leading NA
+ for(i in 1:(width-1)){
+ out[i] <- NA
+ }
+ out
+}
Modified: pkg/GARPFRM/man/EWMA.Rd
===================================================================
--- pkg/GARPFRM/man/EWMA.Rd 2014-03-11 23:04:15 UTC (rev 117)
+++ pkg/GARPFRM/man/EWMA.Rd 2014-03-12 15:40:12 UTC (rev 118)
@@ -32,4 +32,42 @@
correlation by minimizing the sum of squared differences
between the estimated value and realized value.
}
+\examples{
+# data and parameters for EWMA estimate
+data(crsp_weekly)
+R <- largecap_weekly[, 1:2]
+mvR <- largecap_weekly[,1:4]
+lambda <- 0.94
+initialWindow <- 150
+# volatility estimate of univariate data
+lambda <- estimateLambdaVol(R[,1], initialWindow, n=10)
+vol1 <- EWMA(R[,1], lambda=NULL, initialWindow, n=10, "volatility")
+vol1a <- EWMA(R[,1], lambda, initialWindow, n=10, "volatility")
+all.equal(vol1$estimate, vol1a$estimate)
+vol1
+
+# covariance estimate of bivariate data
+lambda <- estimateLambdaCov(R, initialWindow, n=10)
+cov1 <- EWMA(R, lambda=NULL, initialWindow, n=10, "covariance")
+cov1a <- EWMA(R, lambda, initialWindow, n=10, "covariance")
+all.equal(cov1$estimate, cov1a$estimate)
+cov1
+
+# correlation estimate of bivariate data
+lambda <- estimateLambdaCor(R, initialWindow, n=10)
+cor1 <- EWMA(R, lambda=NULL, initialWindow, n=10, "correlation")
+cor1a <- EWMA(R, lambda, initialWindow, n=10, "correlation")
+all.equal(cor1$estimate, cor1a$estimate)
+cor1
+
+# Multivariate EWMA estimate of covariance
+lambda <- 0.94
+cov_mv <- EWMA(mvR, lambda, initialWindow, type="covariance")
+cov_mv
+
+# Multivariate EWMA estimate of correlation
+cor_mv <- EWMA(mvR, lambda, initialWindow, type="correlation")
+cor_mv
+}
+
Modified: pkg/GARPFRM/man/estimateLambdaCor.Rd
===================================================================
--- pkg/GARPFRM/man/estimateLambdaCor.Rd 2014-03-11 23:04:15 UTC (rev 117)
+++ pkg/GARPFRM/man/estimateLambdaCor.Rd 2014-03-12 15:40:12 UTC (rev 118)
@@ -16,4 +16,15 @@
\description{
Estimate lambda for EWMA correlation estimate
}
+\details{
+ The optimal value for lambda is calcualted by minimizing
+ the mean squared error between the estimated correlation
+ and realized correlation.
+}
+\examples{
+data(crsp_weekly)
+R <- largecap_weekly[, 1:2]
+initialWindow <- 150
+lambda <- estimateLambdaCor(R, initialWindow, n=10)
+}
Modified: pkg/GARPFRM/man/estimateLambdaCov.Rd
===================================================================
--- pkg/GARPFRM/man/estimateLambdaCov.Rd 2014-03-11 23:04:15 UTC (rev 117)
+++ pkg/GARPFRM/man/estimateLambdaCov.Rd 2014-03-12 15:40:12 UTC (rev 118)
@@ -16,4 +16,15 @@
\description{
Estimate lambda for EWMA covariance estimate
}
+\details{
+ The optimal value for lambda is calcualted by minimizing
+ the mean squared error between the estimated covariance
+ and realized covariance.
+}
+\examples{
+data(crsp_weekly)
+R <- largecap_weekly[, 1:2]
+initialWindow <- 150
+lambda <- estimateLambdaCov(R, initialWindow, n=10)
+}
Modified: pkg/GARPFRM/man/estimateLambdaVol.Rd
===================================================================
--- pkg/GARPFRM/man/estimateLambdaVol.Rd 2014-03-11 23:04:15 UTC (rev 117)
+++ pkg/GARPFRM/man/estimateLambdaVol.Rd 2014-03-12 15:40:12 UTC (rev 118)
@@ -16,4 +16,15 @@
\description{
Estimate lambda for EWMA volatility estimate
}
+\details{
+ The optimal value for lambda is calcualted by minimizing
+ the mean squared error between the estimated volatility
+ and realized volatility.
+}
+\examples{
+data(crsp_weekly)
+R <- largecap_weekly[, 1]
+initialWindow <- 150
+lambda <- estimateLambdaVol(R, initialWindow, n=10)
+}
Modified: pkg/GARPFRM/man/getCor.Rd
===================================================================
--- pkg/GARPFRM/man/getCor.Rd 2014-03-11 23:04:15 UTC (rev 117)
+++ pkg/GARPFRM/man/getCor.Rd 2014-03-12 15:40:12 UTC (rev 118)
@@ -11,7 +11,16 @@
assets can be specified by name or index.}
}
\description{
- Extract the correlation of two assets from an \code{EWMA}
- object
+ Extract the correlation of two assets from an
+ \code{mvEWMAcor} object
}
+\examples{
+data(crsp_weekly)
+mvR <- largecap_weekly[,1:4]
+lambda <- 0.94
+initialWindow <- 150
+cor_mv <- EWMA(mvR, lambda, initialWindow, type="correlation")
+# Extract the estimated correlation between ORCL and MSFT
+tail(getCor(cor_mv, assets=c("ORCL", "MSFT")))
+}
Modified: pkg/GARPFRM/man/getCov.Rd
===================================================================
--- pkg/GARPFRM/man/getCov.Rd 2014-03-11 23:04:15 UTC (rev 117)
+++ pkg/GARPFRM/man/getCov.Rd 2014-03-12 15:40:12 UTC (rev 118)
@@ -11,7 +11,16 @@
assets can be specified by name or index.}
}
\description{
- Extract the covariance of two assets from an \code{EWMA}
- object
+ Extract the covariance of two assets from an
+ \code{mvEWMAcov} object
}
+\examples{
+data(crsp_weekly)
+mvR <- largecap_weekly[,1:4]
+lambda <- 0.94
+initialWindow <- 150
+cov_mv <- EWMA(mvR, lambda, initialWindow, type="covariance")
+# Extract the estimated covariance between ORCL and MSFT
+tail(getCov(cov_mv, assets=c("ORCL", "MSFT")))
+}
Modified: pkg/GARPFRM/man/plot.EWMA.Rd
===================================================================
--- pkg/GARPFRM/man/plot.EWMA.Rd 2014-03-11 23:04:15 UTC (rev 117)
+++ pkg/GARPFRM/man/plot.EWMA.Rd 2014-03-12 15:40:12 UTC (rev 118)
@@ -14,10 +14,11 @@
\item{\dots}{additional arguments passed to
\code{plot.xts}}
- \item{assets}{for multivariate EWMA objects, character
- vector or numeric vector of assets to extract from the
- covariance or correlation matrix. The assets can be
- specified by name or index.}
+ \item{assets}{character vector or numeric vector of
+ assets to extract from the covariance or correlation
+ matrix. The assets can be specified by name or index.
+ This argument is only usd for multivariate EWMA estimates
+ of a covariance or correlation matrix.}
\item{legendLoc}{location of legend. If NULL, the legend
will be omitted from the plot}
@@ -25,6 +26,38 @@
\item{main}{main title for the plot}
}
\description{
- Plot method for EWMA objects
+ Plot method for EWMA objects.
}
+\examples{
+# data and parameters for EWMA estimate
+data(crsp_weekly)
+R <- largecap_weekly[, 1:2]
+mvR <- largecap_weekly[,1:4]
+lambda <- 0.94
+initialWindow <- 150
+# volatility estimate of univariate data
+vol1 <- EWMA(R[,1], lambda, initialWindow, type="volatility")
+plot(vol1)
+
+# covariance estimate of bivariate data
+cov1 <- EWMA(R, lambda, initialWindow, type="covariance")
+plot(cov1)
+
+# correlation estimate of bivariate data
+cor1 <- EWMA(R, lambda, initialWindow, type="correlation")
+plot(cor1)
+
+# Multivariate EWMA estimate of covariance
+cov_mv <- EWMA(mvR, lambda, initialWindow, type="covariance")
+# These two are equivalent
+plot(cov_mv, assets=c("ORCL", "MSFT"))
+plot(cov_mv, assets=c(1, 2))
+
+# Multivariate EWMA estimate of correlation
+cor_mv <- EWMA(mvR, lambda, initialWindow, type="correlation")
+# These two are equivalent
+plot(cor_mv, assets=c("ORCL", "EMC"))
+plot(cor_mv, assets=c(1, 4))
+}
+
Modified: pkg/GARPFRM/man/realizedCor.Rd
===================================================================
--- pkg/GARPFRM/man/realizedCor.Rd 2014-03-11 23:04:15 UTC (rev 117)
+++ pkg/GARPFRM/man/realizedCor.Rd 2014-03-12 15:40:12 UTC (rev 118)
@@ -13,4 +13,15 @@
\description{
Calculate realized correlation
}
+\details{
+ Realized correlation is calculated using the previous
+ \code{n} periods.
+}
+\examples{
+data(crsp_weekly)
+R <- largecap_weekly[, 1:2]
+# Calculate realized correlation
+realizedCorrelation <- realizedCor(R, 10)
+head(realizedCorrelation)
+}
Modified: pkg/GARPFRM/man/realizedCov.Rd
===================================================================
--- pkg/GARPFRM/man/realizedCov.Rd 2014-03-11 23:04:15 UTC (rev 117)
+++ pkg/GARPFRM/man/realizedCov.Rd 2014-03-12 15:40:12 UTC (rev 118)
@@ -13,4 +13,15 @@
\description{
Calculate realized covariance
}
+\details{
+ Realized covariance is calculated using the previous
+ \code{n} periods.
+}
+\examples{
+data(crsp_weekly)
+R <- largecap_weekly[, 1:2]
+# Calculate realized covariance
+realizedCovariance <- realizedCov(R, 10)
+head(realizedCovariance)
+}
Modified: pkg/GARPFRM/man/realizedVol.Rd
===================================================================
--- pkg/GARPFRM/man/realizedVol.Rd 2014-03-11 23:04:15 UTC (rev 117)
+++ pkg/GARPFRM/man/realizedVol.Rd 2014-03-12 15:40:12 UTC (rev 118)
@@ -13,4 +13,15 @@
\description{
Calculate realized volatility
}
+\details{
+ Realized volatility is defined as the standard deviation
+ of using the previous \code{n} periods.
+}
+\examples{
+data(crsp_weekly)
+R <- largecap_weekly[, 1]
+# Calculate realized volatility
+realizedVolatility <- realizedVol(R[,1], 10)
+head(realizedVolatility)
+}
Added: pkg/GARPFRM/sandbox/loop_rollapply_benchmark.R
===================================================================
--- pkg/GARPFRM/sandbox/loop_rollapply_benchmark.R (rev 0)
+++ pkg/GARPFRM/sandbox/loop_rollapply_benchmark.R 2014-03-12 15:40:12 UTC (rev 118)
@@ -0,0 +1,49 @@
+library(xts)
+
+# benchmark a rolling covariance using a for loop vs rollapply
+
+foo <- function(x, width, na.pad=FALSE){
+ n <- nrow(R)
+ out <- vector("numeric", n-width+1)
+ for(i in width:n){
+ #cat("from = ", i-width+1, "\n")
+ #cat("to = ", i, "\n")
+ tmpR <- R[(i-width+1):i,]
+ out[(i-width+1)] <- cov(tmpR[,1], tmpR[,2])
+ }
+ if(na.pad) {
+ vecNA <- rep(NA, width-1)
+ out <- c(vecNA, out)
+ }
+ out
+}
+
+rollfoo <- function(x, width){
+ rollapply(R, width=width, FUN=function(x) cov(x=x[,1], y=x[,2]), by.column=FALSE)
+}
+
+# set.seed(123)
+# x <- rnorm(15, 0, 0.1)
+# y <- rnorm(15, 0, 0.15)
+# R <- cbind(x, y)
+#
+# cov(R[1:5,])[1,2]
+# cov(R[2:6,])[1,2]
+# cov(R[3:7,])[1,2]
+# cov(R[4:8,])[1,2]
+
+set.seed(123)
+x <- rnorm(1e4, 0, 0.1)
+y <- rnorm(1e4, 0, 0.15)
+R <- cbind(x, y)
+
+cov(R[1:5,])[1,2]
+cov(R[2:6,])[1,2]
+cov(R[3:7,])[1,2]
+cov(R[4:8,])[1,2]
+cov(R[5:9,])[1,2]
+foo(R, 5)[1:5]
+
+all.equal(foo(R, 5), rollfoo(R, 5))
+
+rbenchmark::benchmark(foo(R, 5), rollfoo(R, 5))
Modified: pkg/GARPFRM/sandbox/ross_EWMA.R
===================================================================
--- pkg/GARPFRM/sandbox/ross_EWMA.R 2014-03-11 23:04:15 UTC (rev 117)
+++ pkg/GARPFRM/sandbox/ross_EWMA.R 2014-03-12 15:40:12 UTC (rev 118)
@@ -7,7 +7,7 @@
R <- largecap_weekly[, 1:2]
mvR <- largecap_weekly[,1:4]
lambda <- 0.94
-initialWindow <- 15
+initialWindow <- 150
# volatility estimate of univariate data
lambda <- estimateLambdaVol(R[,1], initialWindow, n=10)
@@ -18,48 +18,53 @@
plot(vol1)
# Calculate realized volatility
-# realizedVol(R, 10)
+realizedVolatility <- realizedVol(R[,1], 10)
# covariance estimate of bivariate data
-lambda <- 0.94
-cov1 <- EWMA(R, lambda, initialWindow, type="covariance")
-cov1a <- uvEWMAcov(R, lambda, initialWindow)
-all.equal(cov1$estimate, cov1a)
+lambda <- estimateLambdaCov(R, initialWindow, n=10)
+cov1 <- EWMA(R, lambda=NULL, initialWindow, n=10, "covariance")
+cov1a <- EWMA(R, lambda, initialWindow, n=10, "covariance")
+all.equal(cov1$estimate, cov1a$estimate)
+cov1
+plot(cov1)
# Calculate realized covariance
-# realizedCov(R, 10)
+realizedCov <- realizedCov(R, 10)
# correlation estimate of bivariate data
-cor1 <- EWMA(R, lambda, initialWindow, type="correlation")
-cor1a <- uvEWMAcor(R, lambda, initialWindow)
-all.equal(cor1$estimate, cor1a)
+lambda <- estimateLambdaCor(R, initialWindow, n=10)
+cor1 <- EWMA(R, lambda=NULL, initialWindow, n=10, "correlation")
+cor1a <- EWMA(R, lambda, initialWindow, n=10, "correlation")
+all.equal(cor1$estimate, cor1a$estimate)
+cor1
+plot(cor1)
-# Calcualte realized correlation
-# realizedCor(R, 10)
+# Calculate realized correlation
+realizedCorrelation <- realizedCor(R, 10)
+# Multivariate EWMA estimate of covariance
+lambda <- 0.94
cov_mv <- EWMA(mvR, lambda, initialWindow, type="covariance")
cov_mv
-tail(getCov(cov_mv, assets=c(1,2)))
-plot(cov_mv)
+# Extract the estimated covariance between ORCL and MSFT
+tail(getCov(cov_mv, assets=c("ORCL", "MSFT")))
+# These two are equivalent
+plot(cov_mv, assets=c("ORCL", "MSFT"))
+plot(cov_mv, assets=c(1, 2))
+
+# Multivariate EWMA estimate of correlation
cor_mv <- EWMA(mvR, lambda, initialWindow, type="correlation")
-tail(getCor(cor_mv, assets=c(1,2)))
+
+# Extract the estimated correlation between ORCL and EMC
+tail(getCor(cor_mv, assets=c("ORCL", "EMC")))
cor_mv
-plot(cor_mv)
-# cor.fun <- function(x){
-# cor(x)[1,2]
-# }
-#
-# cov.fun <- function(x){
-# print(cov(x))
-# cov(x)[1,2]
-# }
-#
-# cor.fun(R)
-# cor(R)
-#
-# cov(R)[1,2]
-#
-# rollapply(data=R, x=R[,1], y=R[,2], width=10, FUN=cov)
-# rollapply(data=R, width=10, FUN=cov.fun)
+# These two are equivalent
+plot(cor_mv, assets=c("ORCL", "EMC"))
+plot(cor_mv, assets=c(1, 4))
+
+
+
+
+
Modified: pkg/GARPFRM/vignettes/EstimatingVolatilitiesCorrelation.Rnw
===================================================================
--- pkg/GARPFRM/vignettes/EstimatingVolatilitiesCorrelation.Rnw 2014-03-11 23:04:15 UTC (rev 117)
+++ pkg/GARPFRM/vignettes/EstimatingVolatilitiesCorrelation.Rnw 2014-03-12 15:40:12 UTC (rev 118)
@@ -103,8 +103,8 @@
<<>>=
lambda <- 0.94
initialWindow <- 15
-volEst <- volEWMA(R, lambda, initialWindow)
-tail(volEst)
+volEst <- EWMA(R, lambda, initialWindow, type="volatility")
+volEst
@
An important point to note is that we are using weekly returns to estimate weekly volatility while the lambda value used in the RiskMetrics database is for daily volatility. A data driven approach for selecting a value for $\lambda$ is to determine the $\lambda$ that minimizes the sum of squared errors between the realized volatility and the estimated volatility from the EWMA model.
@@ -116,35 +116,51 @@
Now we plot the estimated volatility from the EWMA model and the realized volatility.
<<>>=
-plot(vol, main="EWMA Volatility Estimate vs. Realized Volatility")
-lines(volEst, col="red")
+plot(volEst, main="EWMA Volatility Estimate vs. Realized Volatility")
+lines(vol, col="red")
@
-The \code{estimateLambda} function estimates the value for $\lambda$ by minimizing the sum of squared errors between the realized volatility and the EWMA model estimated volatility.
+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.
<<>>=
# Estimate lambda
# Use initialWindow = 15 for the EWMA volatility estimate and
# n = 5 to calculate the realized volatility
-lambda <- estimateLambda(R, initialWindow, n=5)
+lambda <- estimateLambdaVol(R, initialWindow, n=5)
lambda
-volEst2 <- volEWMA(R, lambda, initialWindow)
+volEst2 <- EWMA(R, lambda, initialWindow, type="volatility")
+volEst2
@
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/uwgarp -r 118
More information about the Uwgarp-commits
mailing list