[Uwgarp-commits] r127 - in pkg/GARPFRM: R sandbox
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Sun Mar 23 05:13:36 CET 2014
Author: rossbennett34
Date: 2014-03-23 05:13:33 +0100 (Sun, 23 Mar 2014)
New Revision: 127
Added:
pkg/GARPFRM/R/boot.R
pkg/GARPFRM/R/utils.R
pkg/GARPFRM/sandbox/test_boot.R
Log:
Adding functions for bootstrap
Added: pkg/GARPFRM/R/boot.R
===================================================================
--- pkg/GARPFRM/R/boot.R (rev 0)
+++ pkg/GARPFRM/R/boot.R 2014-03-23 04:13:33 UTC (rev 127)
@@ -0,0 +1,237 @@
+.bootFUN <- function(R, FUN="mean", ..., replications=1000, parallel=FALSE){
+ # R should be a univariate xts object
+
+ fun <- match.fun(FUN)
+ if(!is.function(fun)) stop("FUN could not be matched")
+
+ if(is.function(fun)){
+ .formals <- formals(fun)
+ # add the dots
+ .formals <- modify.args(formals=.formals, ...=..., dots=TRUE)
+ .formals$... <- NULL
+ }
+ # print(.formals)
+
+ replications <- as.integer(replications)
+ n <- nrow(R)
+ out <- vector("numeric", replications)
+
+ if(parallel){
+ stopifnot("package:foreach" %in% search() || require("foreach",quietly = TRUE))
+ out <- foreach(i=1:replications, .inorder=FALSE, .combine=c, .errorhandling='remove') %dopar% {
+ tmpR <- R[sample.int(n, replace=TRUE)]
+ # match the resampled data to R or x in .formals
+ if("R" %in% names(.formals)){
+ .formals <- modify.args(formals=.formals, arglist=NULL, R=tmpR, dots=TRUE)
+ } else if("x" %in% names(.formals)){
+ .formals <- modify.args(formals=.formals, arglist=NULL, x=tmpR, dots=TRUE)
+ }
+ do.call(fun, .formals)
+ }
+ } else {
+ for(i in 1:replications){
+ # sampled data
+ tmpR <- R[sample.int(n, replace=TRUE),]
+ # match the resampled data to R or x in .formals
+ if("R" %in% names(.formals)){
+ .formals <- modify.args(formals=.formals, arglist=NULL, R=tmpR, dots=TRUE)
+ } else if("x" %in% names(.formals)){
+ .formals <- modify.args(formals=.formals, arglist=NULL, x=tmpR, dots=TRUE)
+ }
+ # call the function
+ tmp <- try(do.call(fun, .formals), silent=TRUE)
+
+ # if try-error, stop the function call, else insert the resampled statistic
+ # to the output vector
+ if(inherits(tmp, "try-error")){
+ stop("FUN could not be evaluated")
+ } else {
+ out[i] <- tmp
+ }
+ }
+ }
+ # compute the expected value of the statistic on resampled data
+ mean(out)
+}
+
+.bootMean <- function(R, ..., replications=1000, parallel=FALSE){
+ # R should be a univariate xts object
+ .bootFUN(R=R, FUN="mean", ...=..., replications=replications, parallel=parallel)
+}
+
+bootMean <- function(R, ..., replications=1000, parallel=FALSE){
+ if(ncol(R) == 1){
+ tmp <- .bootMean(R=R, ...=..., replications=replications, parallel=parallel)
+ } else {
+ tmp <- vector("numeric", ncol(R))
+ for(i in 1:ncol(R)){
+ tmp[i] <- .bootMean(R=R[,i], ...=..., replications=replications, parallel=parallel)
+ }
+ }
+ out <- matrix(tmp, nrow=1, ncol=ncol(R))
+ rownames(out) <- "mean"
+ colnames(out) <- colnames(R)
+ return(out)
+}
+
+.bootSD <- function(R, ..., replications=1000, parallel=FALSE){
+ # R should be a univariate xts object
+ .bootFUN(R=R, FUN="sd", ...=..., replications=replications, parallel=parallel)
+}
+
+bootSD <- function(R, ..., replications=1000, parallel=FALSE){
+ if(ncol(R) == 1){
+ tmp <- .bootSD(R=R, ...=..., replications=replications, parallel=parallel)
+ } else {
+ tmp <- vector("numeric", ncol(R))
+ for(i in 1:ncol(R)){
+ tmp[i] <- .bootSD(R=R[,i], ...=..., replications=replications, parallel=parallel)
+ }
+ }
+ out <- matrix(tmp, nrow=1, ncol=ncol(R))
+ rownames(out) <- "sd"
+ colnames(out) <- colnames(R)
+ return(out)
+}
+
+.bootStdDev <- function(R, ..., replications=1000, parallel=FALSE){
+ # R should be a univariate xts object
+ .bootFUN(R=R, FUN="StdDev", ...=..., replications=replications, parallel=parallel)
+}
+
+bootStdDev <- function(R, ..., replications=1000, parallel=FALSE){
+ if(ncol(R) == 1){
+ tmp <- .bootStdDev(R=R, ...=..., replications=replications, parallel=parallel)
+ } else {
+ tmp <- vector("numeric", ncol(R))
+ for(i in 1:ncol(R)){
+ tmp[i] <- .bootStdDev(R=R[,i], ...=..., replications=replications, parallel=parallel)
+ }
+ }
+ out <- matrix(tmp, nrow=1, ncol=ncol(R))
+ rownames(out) <- "StdDev"
+ colnames(out) <- colnames(R)
+ return(out)
+}
+
+tmpCor <- function(R){
+ # R should be a bivariate xts object
+ cor(x=R[,1], y=R[,2])
+}
+
+.bootCor <- function(R, ..., replications=1000, parallel=FALSE){
+ # R should be a bivariate xts object
+ .bootFUN(R=R[,1:2], FUN="tmpCor", ...=..., replications=replications, parallel=parallel)
+}
+
+bootCor <- function(R, ..., replications=1000, parallel=FALSE){
+ cnames <- colnames(R)
+ if(ncol(R) == 2){
+ tmp <- .bootCor(R=R, ...=..., replications=replications, parallel=parallel)
+ out_names <- paste(cnames[1], cnames[2], sep=".")
+ num_col <- 1
+ } else {
+ tmp <- vector("numeric", choose(ncol(R), 2))
+ out_names <- vector("numeric", length(tmp))
+ num_col <- length(tmp)
+ k <- 1
+ for(i in 1:(ncol(R)-1)){
+ for(j in (i+1):ncol(R)){
+ tmp[k] <- .bootCor(R=cbind(R[,i], R[,j]), ...=..., replications=replications, parallel=parallel)
+ out_names[k] <- paste(cnames[i], cnames[j], sep=".")
+ k <- k + 1
+ }
+ }
+ }
+ # out <- matrix(tmp, nrow=1, ncol=ncol(R))
+ out <- matrix(tmp, nrow=1, ncol=num_col)
+ rownames(out) <- "cor"
+ colnames(out) <- out_names
+ return(out)
+}
+
+tmpCov <- function(R){
+ # R should be a bivariate xts object
+ cov(x=R[,1], y=R[,2])
+}
+
+.bootCov <- function(R, ..., replications=1000, parallel=FALSE){
+ # R should be a bivariate xts object
+ .bootFUN(R=R[,1:2], FUN="tmpCov", ...=..., replications=replications, parallel=parallel)
+}
+
+bootCov <- function(R, ..., replications=1000, parallel=FALSE){
+ cnames <- colnames(R)
+ if(ncol(R) == 2){
+ tmp <- .bootCov(R=R, ...=..., replications=replications, parallel=parallel)
+ out_names <- paste(cnames[1], cnames[2], sep=".")
+ num_col <- 1
+ } else {
+ tmp <- vector("numeric", choose(ncol(R), 2))
+ out_names <- vector("numeric", length(tmp))
+ num_col <- length(tmp)
+ k <- 1
+ for(i in 1:(ncol(R)-1)){
+ for(j in (i+1):ncol(R)){
+ tmp[k] <- .bootCov(R=cbind(R[,i], R[,j]), ...=..., replications=replications, parallel=parallel)
+ out_names[k] <- paste(cnames[i], cnames[j], sep=".")
+ k <- k + 1
+ }
+ }
+ }
+ # out <- matrix(tmp, nrow=1, ncol=ncol(R))
+ out <- matrix(tmp, nrow=1, ncol=num_col)
+ rownames(out) <- "cov"
+ colnames(out) <- out_names
+ return(out)
+}
+
+.bootVaR <- function(R, ..., replications=1000, parallel=FALSE){
+ # R should be a univariate xts object
+ .bootFUN(R=R[,1], FUN="VaR", ...=..., replications=replications, parallel=parallel)
+}
+
+bootVaR <- function(R, ..., replications=1000, parallel=FALSE){
+ if(ncol(R) == 1){
+ tmp <- .bootVaR(R=R, ...=..., replications=replications, parallel=parallel)
+ } else {
+ tmp <- vector("numeric", ncol(R))
+ for(i in 1:ncol(R)){
+ tmp[i] <- .bootVaR(R=R[,i], ...=..., replications=replications, parallel=parallel)
+ }
+ }
+ out <- matrix(tmp, nrow=1, ncol=ncol(R))
+ rownames(out) <- "VaR"
+ colnames(out) <- colnames(R)
+ return(out)
+}
+
+.bootES <- function(R, ..., replications=1000, parallel=FALSE){
+ # R should be a univariate xts object
+ .bootFUN(R=R, FUN="ES", ...=..., replications=replications, parallel=parallel)
+}
+
+bootES <- function(R, ..., replications=1000, parallel=FALSE){
+ if(ncol(R) == 1){
+ tmp <- .bootES(R=R, ...=..., replications=replications, parallel=parallel)
+ } else {
+ tmp <- vector("numeric", ncol(R))
+ for(i in 1:ncol(R)){
+ tmp[i] <- .bootES(R=R[,i], ...=..., replications=replications, parallel=parallel)
+ }
+ }
+ out <- matrix(tmp, nrow=1, ncol=ncol(R))
+ rownames(out) <- "ES"
+ colnames(out) <- colnames(R)
+ return(out)
+}
+
+# .bootMean <- function(R, replications=1000){
+# replications <- as.integer(replications)
+# n <- length(R)
+# out <- vector("numeric", replications)
+# for(i in 1:replications){
+# out[i] <- mean(R[sample.int(n, replace=TRUE)])
+# }
+# mean(out)
+# }
Added: pkg/GARPFRM/R/utils.R
===================================================================
--- pkg/GARPFRM/R/utils.R (rev 0)
+++ pkg/GARPFRM/R/utils.R 2014-03-23 04:13:33 UTC (rev 127)
@@ -0,0 +1,38 @@
+
+modify.args <- function(formals, arglist, ..., dots=FALSE){
+ # modify.args function from quantstrat
+
+ # avoid evaluating '...' to make things faster
+ dots.names <- eval(substitute(alist(...)))
+
+ if(missing(arglist))
+ arglist <- NULL
+ arglist <- c(arglist, dots.names)
+
+ # see 'S Programming' p. 67 for this matching
+
+ # nothing to do if arglist is empty; return formals
+ if(!length(arglist))
+ return(formals)
+
+ argnames <- names(arglist)
+ if(!is.list(arglist) && !is.null(argnames) && !any(argnames == ""))
+ stop("'arglist' must be a *named* list, with no names == \"\"")
+
+ .formals <- formals
+ onames <- names(.formals)
+
+ pm <- pmatch(argnames, onames, nomatch = 0L)
+ #if(any(pm == 0L))
+ # message(paste("some arguments stored for", fun, "do not match"))
+ names(arglist[pm > 0L]) <- onames[pm]
+ .formals[pm] <- arglist[pm > 0L]
+
+ # include all elements from arglist if function formals contain '...'
+ if(dots && !is.null(.formals$...)) {
+ dotnames <- names(arglist[pm == 0L])
+ .formals[dotnames] <- arglist[dotnames]
+ #.formals$... <- NULL # should we assume we matched them all?
+ }
+ .formals
+}
Added: pkg/GARPFRM/sandbox/test_boot.R
===================================================================
--- pkg/GARPFRM/sandbox/test_boot.R (rev 0)
+++ pkg/GARPFRM/sandbox/test_boot.R 2014-03-23 04:13:33 UTC (rev 127)
@@ -0,0 +1,65 @@
+# bootstrap
+
+data(crsp_weekly)
+R <- largecap_weekly[,1:4]
+R1 <- R[1:100,1]
+
+set.seed(123)
+.bootFUN(R1, FUN="mean", replications=10000, parallel=FALSE)
+set.seed(123)
+.bootFUN(R1, FUN="mean", replications=10000, parallel=TRUE)
+
+# bootstrap various statistics
+# mean
+bootMean(R[,1])
+bootMean(R)
+
+# sd
+bootSD(R[,1])
+bootSD(R)
+
+# StdDev
+bootStdDev(R[,1])
+bootStdDev(R)
+
+# simpleVolatility
+
+
+# cor
+bootCor(R[,1:2])
+bootCor(R)
+
+# cov
+bootCov(R[,1:2])
+bootCov(R)
+
+# VaR
+bootVaR(R[,1], p=0.9, method="historical")
+bootVaR(R[,1], p=0.9, method="gaussian")
+bootVaR(R, p=0.9, method="historical")
+
+# ES
+bootES(R[,1], p=0.9, method="historical")
+bootES(R, p=0.9, method="historical")
+
+
+# maybe...
+# use bootstrapped returns to imply the prices
+
+# foo1 <- function(x){
+# # Use sample.int and subset
+# x[sample.int(length(x), replace=TRUE)]
+# }
+
+# foo2 <- function(x){
+# # sample directly from the returns
+# sample(x, length(x), replace=TRUE)
+# }
+
+# which is faster?
+
+# set.seed(123)
+# x <- rnorm(1e6)
+# rbenchmark::benchmark(foo1(R[,1]), foo2(R[,1]), replications=1e5)
+
+# foo1 is slightly faster
\ No newline at end of file
More information about the Uwgarp-commits
mailing list