[Returnanalytics-commits] r4027 - in pkg/FactorAnalytics: . R man

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Sat Aug 6 20:03:53 CEST 2016


Author: pragnya
Date: 2016-08-06 20:03:52 +0200 (Sat, 06 Aug 2016)
New Revision: 4027

Modified:
   pkg/FactorAnalytics/DESCRIPTION
   pkg/FactorAnalytics/R/fmEsDecomp.R
   pkg/FactorAnalytics/R/fmVaRDecomp.R
   pkg/FactorAnalytics/man/fmEsDecomp.Rd
   pkg/FactorAnalytics/man/fmVaRDecomp.Rd
Log:
Add option, examples for user specified factor cov and Normal VaR and ES

Modified: pkg/FactorAnalytics/DESCRIPTION
===================================================================
--- pkg/FactorAnalytics/DESCRIPTION	2016-08-05 15:21:42 UTC (rev 4026)
+++ pkg/FactorAnalytics/DESCRIPTION	2016-08-06 18:03:52 UTC (rev 4027)
@@ -1,8 +1,8 @@
 Package: factorAnalytics
 Type: Package
 Title: Factor Analytics
-Version: 2.0.36
-Date: 2016-08-04
+Version: 2.0.37
+Date: 2016-08-06
 Author: Eric Zivot, Sangeetha Srinivasan and Yi-An Chen
 Maintainer: Sangeetha Srinivasan <sangee at uw.edu>
 Description: Linear factor model fitting for asset returns (three major types-

Modified: pkg/FactorAnalytics/R/fmEsDecomp.R
===================================================================
--- pkg/FactorAnalytics/R/fmEsDecomp.R	2016-08-05 15:21:42 UTC (rev 4026)
+++ pkg/FactorAnalytics/R/fmEsDecomp.R	2016-08-06 18:03:52 UTC (rev 4027)
@@ -19,22 +19,22 @@
 #' being less than or equal to \code{VaR.fm}. This is estimated as a sample 
 #' average of the observations in that data window. 
 #' 
+#' Refer to Eric Zivot's slides (referenced) for formulas pertaining to the 
+#' calculation of Normal ES (adapted from a portfolio context to factor models).
+#' 
 #' @param object fit object of class \code{tsfm}, \code{sfm} or \code{ffm}.
+#' @param factor.cov optional user specified factor covariance matrix with 
+#' named columns; defaults to the sample covariance matrix.
 #' @param p confidence level for calculation. Default is 0.95.
 #' @param type one of "np" (non-parametric) or "normal" for calculating VaR. 
 #' Default is "np".
-#' @param use an optional character string giving a method for computing factor
-#' covariances in the presence of missing values. This must be (an 
-#' abbreviation of) one of the strings "everything", "all.obs", 
-#' "complete.obs", "na.or.complete", or "pairwise.complete.obs". Default is 
-#' "pairwise.complete.obs".
+#' @param use method for computing covariances in the presence of missing 
+#' values; one of "everything", "all.obs", "complete.obs", "na.or.complete", or 
+#' "pairwise.complete.obs". Default is "pairwise.complete.obs".
 #' @param ... other optional arguments passed to \code{\link[stats]{quantile}}.
 #' 
 #' @return A list containing 
 #' \item{ES.fm}{length-N vector of factor model ES of N-asset returns.}
-#' \item{n.exceed}{length-N vector of number of observations beyond VaR for 
-#' each asset.}
-#' \item{idx.exceed}{list of numeric vector of index values of exceedances.}
 #' \item{mES}{N x (K+1) matrix of marginal contributions to VaR.}
 #' \item{cES}{N x (K+1) matrix of component contributions to VaR.}
 #' \item{pcES}{N x (K+1) matrix of percentage component contributions to VaR.}
@@ -43,6 +43,9 @@
 #' @author Eric Zviot, Sangeetha Srinivasan and Yi-An Chen
 #' 
 #' @references 
+#' Eric Zivot's slides from CFRM 546: Estimating risk measures: Portfolio of 
+#' Assets, April 28, 2015.
+#' 
 #' Epperlein, E., & Smillie, A. (2006). Portfolio risk analysis Cracking VAR 
 #' with kernels. RISK-LONDON-RISK MAGAZINE LIMITED-, 19(8), 70.
 #' 
@@ -63,7 +66,7 @@
 #' \code{\link{fmVaRDecomp}} for factor model VaR decomposition.
 #' 
 #' @examples
-#' # Time Series Factor Model
+#' #' # Time Series Factor Model
 #' data(managers)
 #' fit.macro <- fitTsfm(asset.names=colnames(managers[,(1:6)]),
 #'                      factor.names=colnames(managers[,(7:8)]), data=managers)
@@ -74,9 +77,17 @@
 #' # Statistical Factor Model
 #' data(StockReturns)
 #' sfm.pca.fit <- fitSfm(r.M, k=2)
-#' ES.decomp <- fmEsDecomp(sfm.pca.fit)
+#' ES.decomp <- fmEsDecomp(sfm.pca.fit, type="normal")
 #' ES.decomp$cES
 #' 
+#' # Fundamental Factor Model
+#' data(Stock.df)
+#' exposure.vars <- c("BOOK2MARKET", "LOG.MARKETCAP")
+#' fit <- fitFfm(data=stock, asset.var="TICKER", ret.var="RETURN", 
+#'               date.var="DATE", exposure.vars=exposure.vars)
+#' ES.decomp <- fmEsDecomp(fit, type="normal")
+#' head(ES.decomp$cES)
+#' 
 #' @export
 
 fmEsDecomp <- function(object, ...){
@@ -91,9 +102,8 @@
 #' @method fmEsDecomp tsfm
 #' @export
 
-fmEsDecomp.tsfm <- function(object, p=0.95, type=c("np","normal"),  
+fmEsDecomp.tsfm <- function(object, factor.cov, p=0.95, type=c("np","normal"),  
                             use="pairwise.complete.obs", ...) {
-  
   # set default for type
   type = type[1]
   
@@ -112,14 +122,36 @@
   resid.xts <- as.xts(t(t(residuals(object))/object$resid.sd))
   time(resid.xts) <- as.Date(time(resid.xts))
   
+  if (type=="normal") {
+    # get cov(F): K x K
+    if (missing(factor.cov)) {
+      factor.cov = cov(as.matrix(factors.xts), use=use, ...)
+    } else {
+      if (!identical(dim(factor.cov), as.integer(c(ncol(factor), ncol(factor))))) {
+        stop("Dimensions of user specified factor covariance matrix are not
+           compatible with the number of factors in the fitTsfm object")
+      }
+    }
+    # get cov(F.star): (K+1) x (K+1)
+    K <- ncol(object$beta)
+    factor.star.cov <- diag(K+1)
+    factor.star.cov[1:K, 1:K] <- factor.cov
+    colnames(factor.star.cov) <- c(colnames(factor.cov),"residuals")
+    rownames(factor.star.cov) <- c(colnames(factor.cov),"residuals")
+    # factor expected returns
+    MU <- c(colMeans(factors.xts, na.rm=TRUE), 0)
+    names(MU) <- colnames(beta.star)
+    # SIGMA*Beta to compute normal mES
+    SIGB <-  beta.star %*% factor.star.cov
+  }
+  
   # initialize lists and matrices
   N <- length(object$asset.names)
   K <- length(object$factor.names)
   VaR.fm <- rep(NA, N)
   ES.fm <- rep(NA, N)
   idx.exceed <- list()
-  n.exceed <- rep(NA, N)
-  names(VaR.fm) = names(ES.fm) = names(n.exceed) = object$asset.names
+  names(VaR.fm) = names(ES.fm) = object$asset.names
   mES <- matrix(NA, N, K+1)
   cES <- matrix(NA, N, K+1)
   pcES <- matrix(NA, N, K+1)
@@ -129,29 +161,31 @@
   for (i in object$asset.names) {
     # return data for asset i
     R.xts <- object$data[,i]
-    # get VaR for asset i
-    if (type=="np") {
+    
+    if (type=="np") { 
+      # get VaR for asset i
       VaR.fm[i] <- quantile(R.xts, probs=1-p, na.rm=TRUE, ...)
-    } else {
-      VaR.fm[i] <- mean(R.xts, na.rm=TRUE) + sd(R.xts, na.rm=TRUE)*qnorm(1-p)
+      # index of VaR exceedances
+      idx.exceed[[i]] <- which(R.xts <= VaR.fm[i])
+      # compute ES as expected value of asset return, such that the given asset 
+      # return is less than or equal to its value-at-risk (VaR)
+      ES.fm[i] <- mean(R.xts[idx.exceed[[i]]], na.rm =TRUE)
+      # get F.star data object
+      factor.star <- merge(factors.xts, resid.xts[,i])
+      colnames(factor.star)[dim(factor.star)[2]] <- "residual"
+      # compute marginal ES as expected value of factor returns, when the asset's 
+      # return is less than or equal to its value-at-risk (VaR)
+      mES[i,] <- colMeans(factor.star[idx.exceed[[i]],], na.rm =TRUE)
+      
+    } else if (type=="normal") {
+      # extract vector of factor model loadings for asset i
+      beta.i <- beta.star[i,,drop=F]
+      # compute ES
+      ES.fm[i] <- beta.star[i,] %*% MU + sqrt(beta.i %*% factor.star.cov %*% t(beta.i))*dnorm(qnorm(1-p))/(1-p) 
+      # compute marginal ES
+      mES[i,] <- t(MU) + SIGB[i,]/sd(R.xts, na.rm=TRUE) * dnorm(qnorm(1-p))/(1-p)
     }
-    # index of VaR exceedances
-    idx.exceed[[i]] <- which(R.xts <= VaR.fm[i])
-    # number of VaR exceedances
-    n.exceed[i] <- length(idx.exceed[[i]])
     
-    # compute ES as expected value of asset return, such that the given asset 
-    # return is less than or equal to its value-at-risk (VaR)
-    ES.fm[i] <- mean(R.xts[idx.exceed[[i]]], na.rm =TRUE)
-    
-    # get F.star data object
-    factor.star <- merge(factors.xts, resid.xts[,i])
-    colnames(factor.star)[dim(factor.star)[2]] <- "residual"
-    
-    # compute marginal ES as expected value of factor returns, when the asset's 
-    # return is less than or equal to its value-at-risk (VaR)
-    mES[i,] <- colMeans(factor.star[idx.exceed[[i]],], na.rm =TRUE)
-    
     # correction factor to ensure that sum(cES) = asset ES
     cf <- as.numeric( ES.fm[i] / sum(mES[i,]*beta.star[i,], na.rm=TRUE) )
     
@@ -162,8 +196,7 @@
     pcES[i,] <- 100* cES[i,] / ES.fm[i]
   }
   
-  fm.ES.decomp <- list(ES.fm=ES.fm, n.exceed=n.exceed, idx.exceed=idx.exceed, 
-                       mES=mES, cES=cES, pcES=pcES)
+  fm.ES.decomp <- list(ES.fm=ES.fm, mES=mES, cES=cES, pcES=pcES)
   
   return(fm.ES.decomp)
 }
@@ -172,9 +205,8 @@
 #' @method fmEsDecomp sfm
 #' @export
 
-fmEsDecomp.sfm <- function(object, p=0.95, type=c("np","normal"),  
+fmEsDecomp.sfm <- function(object, factor.cov, p=0.95, type=c("np","normal"),  
                            use="pairwise.complete.obs", ...) {
-  
   # set default for type
   type = type[1]
   
@@ -193,46 +225,171 @@
   resid.xts <- as.xts(t(t(residuals(object))/object$resid.sd))
   time(resid.xts) <- as.Date(time(resid.xts))
   
+  if (type=="normal") {
+    # get cov(F): K x K
+    if (missing(factor.cov)) {
+      factor.cov = cov(as.matrix(factors.xts), use=use, ...) 
+    } else {
+      if (!identical(dim(factor.cov), as.integer(c(object$k, object$k)))) {
+        stop("Dimensions of user specified factor covariance matrix are not 
+             compatible with the number of factors in the fitSfm object")
+      }
+    }
+    # get cov(F.star): (K+1) x (K+1)
+    K <- object$k
+    factor.star.cov <- diag(K+1)
+    factor.star.cov[1:K, 1:K] <- factor.cov
+    colnames(factor.star.cov) <- c(colnames(factor.cov),"residuals")
+    rownames(factor.star.cov) <- c(colnames(factor.cov),"residuals")
+    # factor expected returns
+    MU <- c(colMeans(factors.xts, na.rm=TRUE), 0)
+    # SIGMA*Beta to compute normal mVaR
+    SIGB <- beta.star %*% factor.star.cov
+  }
+  
   # initialize lists and matrices
   N <- length(object$asset.names)
   K <- object$k
   VaR.fm <- rep(NA, N)
   ES.fm <- rep(NA, N)
   idx.exceed <- list()
-  n.exceed <- rep(NA, N)
-  names(VaR.fm) = names(ES.fm) = names(n.exceed) = object$asset.names
+  names(VaR.fm) = names(ES.fm) = object$asset.names
   mES <- matrix(NA, N, K+1)
   cES <- matrix(NA, N, K+1)
   pcES <- matrix(NA, N, K+1)
   rownames(mES)=rownames(cES)=rownames(pcES)=object$asset.names
-  colnames(mES)=colnames(cES)=colnames(pcES)=c(paste("F",1:K,sep="."),
-                                               "residuals")
+  colnames(mES)=colnames(cES)=colnames(pcES)=c(paste("F",1:K,sep="."),"residuals")
   
   for (i in object$asset.names) {
     # return data for asset i
     R.xts <- object$data[,i]
-    # get VaR for asset i
+    
     if (type=="np") {
+      # get VaR for asset i
       VaR.fm[i] <- quantile(R.xts, probs=1-p, na.rm=TRUE, ...)
-    } else {
-      VaR.fm[i] <- mean(R.xts, na.rm=TRUE) + sd(R.xts, na.rm=TRUE)*qnorm(1-p)
+      # index of VaR exceedances
+      idx.exceed[[i]] <- which(R.xts <= VaR.fm[i])
+      # compute ES as expected value of asset return, such that the given asset 
+      # return is less than or equal to its value-at-risk (VaR)
+      ES.fm[i] <- mean(R.xts[idx.exceed[[i]]], na.rm =TRUE)
+      # get F.star data object
+      factor.star <- merge(factors.xts, resid.xts[,i])
+      colnames(factor.star)[dim(factor.star)[2]] <- "residual"
+      # compute marginal ES as expected value of factor returns, when the asset's 
+      # return is less than or equal to its value-at-risk (VaR)
+      mES[i,] <- colMeans(factor.star[idx.exceed[[i]],], na.rm =TRUE)
+      
+    } else if (type=="normal") {
+      # extract vector of factor model loadings for asset i
+      beta.i <- beta.star[i,,drop=F]
+      # compute ES
+      ES.fm[i] <- beta.star[i,] %*% MU + sqrt(beta.i %*% factor.star.cov %*% t(beta.i))*dnorm(qnorm(1-p))/(1-p) 
+      # compute marginal ES
+      mES[i,] <- t(MU) + SIGB[i,]/sd(R.xts, na.rm=TRUE) * dnorm(qnorm(1-p))/(1-p)
     }
-    # index of VaR exceedances
-    idx.exceed[[i]] <- which(R.xts <= VaR.fm[i])
-    # number of VaR exceedances
-    n.exceed[i] <- length(idx.exceed[[i]])
     
-    # compute ES as expected value of asset return, such that the given asset 
-    # return is less than or equal to its value-at-risk (VaR)
-    ES.fm[i] <- mean(R.xts[idx.exceed[[i]]], na.rm =TRUE)
+    # correction factor to ensure that sum(cES) = asset ES
+    cf <- as.numeric( ES.fm[i] / sum(mES[i,]*beta.star[i,], na.rm=TRUE) )
     
-    # get F.star data object
-    factor.star <- merge(factors.xts, resid.xts[,i])
-    colnames(factor.star)[dim(factor.star)[2]] <- "residual"
+    # compute marginal, component and percentage contributions to ES
+    # each of these have dimensions: N x (K+1)
+    mES[i,] <- cf * mES[i,]
+    cES[i,] <- mES[i,] * beta.star[i,]
+    pcES[i,] <- 100* cES[i,] / ES.fm[i]
+  }
+  
+  fm.ES.decomp <- list(ES.fm=ES.fm, mES=mES, cES=cES, pcES=pcES)
+  
+  return(fm.ES.decomp)
+}
+
+#' @rdname fmEsDecomp
+#' @method fmEsDecomp ffm
+#' @export
+
+fmEsDecomp.ffm <- function(object, factor.cov, p=0.95, type=c("np","normal"),  
+                           use="pairwise.complete.obs", ...) {
+  # set default for type
+  type = type[1]
+  if (!(type %in% c("np","normal"))) {
+    stop("Invalid args: type must be 'np' or 'normal' ")
+  }
+  
+  # get beta.star
+  beta <- object$beta
+  beta[is.na(beta)] <- 0
+  beta.star <- as.matrix(cbind(beta, sqrt(object$resid.var)))
+  colnames(beta.star)[dim(beta.star)[2]] <- "residual"
+  
+  # factor returns and residuals data
+  factors.xts <- object$factor.returns
+  resid.xts <- as.xts(t(t(residuals(object))/sqrt(object$resid.var)))
+  time(resid.xts) <- as.Date(time(resid.xts))
+  
+  if (type=="normal") {
+    # get cov(F): K x K
+    if (missing(factor.cov)) {
+      factor.cov <- object$factor.cov
+    } else {
+      if (!identical(dim(factor.cov), dim(object$factor.cov))) {
+        stop("Dimensions of user specified factor covariance matrix are not 
+             compatible with the number of factors in the fitSfm object")
+      }
+    }
+    # get cov(F.star): (K+1) x (K+1)
+    K <- ncol(object$beta)
+    factor.star.cov <- diag(K+1)
+    factor.star.cov[1:K, 1:K] <- factor.cov
+    colnames(factor.star.cov) <- c(colnames(factor.cov),"residuals")
+    rownames(factor.star.cov) <- c(colnames(factor.cov),"residuals")
+    # factor expected returns
+    MU <- c(colMeans(factors.xts, na.rm=TRUE), 0)
+    # SIGMA*Beta to compute normal mVaR
+    SIGB <-  beta.star %*% factor.star.cov
+  }
+  
+  # initialize lists and matrices
+  N <- length(object$asset.names)
+  K <- length(object$factor.names)
+  VaR.fm <- rep(NA, N)
+  ES.fm <- rep(NA, N)
+  idx.exceed <- list()
+  names(VaR.fm) = names(ES.fm) = object$asset.names
+  mES <- matrix(NA, N, K+1)
+  cES <- matrix(NA, N, K+1)
+  pcES <- matrix(NA, N, K+1)
+  rownames(mES)=rownames(cES)=rownames(pcES)=object$asset.names
+  colnames(mES)=colnames(cES)=colnames(pcES)=c(object$factor.names, "residuals")
+  
+  for (i in object$asset.names) {
+    # return data for asset i
+    subrows <- which(object$data[[object$asset.var]]==i)
+    R.xts <- as.xts(object$data[subrows,object$ret.var], 
+                    as.Date(object$data[subrows,object$date.var]))
     
-    # compute marginal ES as expected value of factor returns, when the asset's 
-    # return is less than or equal to its value-at-risk (VaR)
-    mES[i,] <- colMeans(factor.star[idx.exceed[[i]],], na.rm =TRUE)
+    if (type=="np") {
+      # get VaR for asset i
+      VaR.fm[i] <- quantile(R.xts, probs=1-p, na.rm=TRUE, ...)
+      # index of VaR exceedances
+      idx.exceed[[i]] <- which(R.xts <= VaR.fm[i])
+      # compute ES as expected value of asset return, such that the given asset 
+      # return is less than or equal to its value-at-risk (VaR)
+      ES.fm[i] <- mean(R.xts[idx.exceed[[i]]], na.rm =TRUE)
+      # get F.star data object
+      factor.star <- merge(factors.xts, resid.xts[,i])
+      colnames(factor.star)[dim(factor.star)[2]] <- "residual"
+      # compute marginal ES as expected value of factor returns, when the asset's 
+      # return is less than or equal to its value-at-risk (VaR)
+      mES[i,] <- colMeans(factor.star[idx.exceed[[i]],], na.rm =TRUE)
+      
+    } else if (type=="normal")  {
+      # extract vector of factor model loadings for asset i
+      beta.i <- beta.star[i,,drop=F]
+      # compute ES
+      ES.fm[i] <- beta.star[i,] %*% MU + sqrt(beta.i %*% factor.star.cov %*% t(beta.i))*dnorm(qnorm(1-p))/(1-p) 
+      # compute marginal ES
+      mES[i,] <- t(MU) + SIGB[i,]/sd(R.xts, na.rm=TRUE) * dnorm(qnorm(1-p))/(1-p)
+    }
     
     # correction factor to ensure that sum(cES) = asset ES
     cf <- as.numeric( ES.fm[i] / sum(mES[i,]*beta.star[i,], na.rm=TRUE) )
@@ -244,8 +401,7 @@
     pcES[i,] <- 100* cES[i,] / ES.fm[i]
   }
   
-  fm.ES.decomp <- list(ES.fm=ES.fm, n.exceed=n.exceed, idx.exceed=idx.exceed, 
-                       mES=mES, cES=cES, pcES=pcES)
+  fm.ES.decomp <- list(ES.fm=ES.fm, mES=mES, cES=cES, pcES=pcES)
   
   return(fm.ES.decomp)
 }

Modified: pkg/FactorAnalytics/R/fmVaRDecomp.R
===================================================================
--- pkg/FactorAnalytics/R/fmVaRDecomp.R	2016-08-05 15:21:42 UTC (rev 4026)
+++ pkg/FactorAnalytics/R/fmVaRDecomp.R	2016-08-06 18:03:52 UTC (rev 4027)
@@ -69,7 +69,6 @@
 #' data(managers)
 #' fit.macro <- fitTsfm(asset.names=colnames(managers[,(1:6)]),
 #'                      factor.names=colnames(managers[,(7:8)]), data=managers)
-#'                      
 #' VaR.decomp <- fmVaRDecomp(fit.macro)
 #' # get the component contributions
 #' VaR.decomp$cVaR
@@ -77,7 +76,6 @@
 #' # Statistical Factor Model
 #' data(StockReturns)
 #' sfm.pca.fit <- fitSfm(r.M, k=2)
-#' 
 #' VaR.decomp <- fmVaRDecomp(sfm.pca.fit, type="normal")
 #' VaR.decomp$cVaR
 #' 
@@ -86,7 +84,6 @@
 #' exposure.vars <- c("BOOK2MARKET", "LOG.MARKETCAP")
 #' fit <- fitFfm(data=stock, asset.var="TICKER", ret.var="RETURN", 
 #'               date.var="DATE", exposure.vars=exposure.vars)
-#' 
 #' VaR.decomp <- fmVaRDecomp(fit, type="normal")
 #' VaR.decomp$cVaR
 #' 
@@ -106,7 +103,6 @@
 
 fmVaRDecomp.tsfm <- function(object, factor.cov, p=0.95, type=c("np","normal"), 
                              use="pairwise.complete.obs", ...) {
-  
   # set default for type
   type = type[1]
   if (!(type %in% c("np","normal"))) {
@@ -134,18 +130,15 @@
            compatible with the number of factors in the fitTsfm object")
       }
     }
-    
     # get cov(F.star): (K+1) x (K+1)
     K <- ncol(object$beta)
     factor.star.cov <- diag(K+1)
     factor.star.cov[1:K, 1:K] <- factor.cov
     colnames(factor.star.cov) <- c(colnames(factor.cov),"residuals")
     rownames(factor.star.cov) <- c(colnames(factor.cov),"residuals")
-    
     # factor expected returns
     MU <- c(colMeans(factors.xts, na.rm=TRUE), 0)
     names(MU) <- colnames(beta.star)
-    
     # SIGMA*Beta to compute normal mVaR
     SIGB <-  beta.star %*% factor.star.cov
   }
@@ -166,20 +159,11 @@
   for (i in object$asset.names) {
     # return data for asset i
     R.xts <- object$data[,i]
-    # get VaR for asset i
-    if (type=="np") {
-      VaR.fm[i] <- quantile(R.xts, probs=1-p, na.rm=TRUE, ...)
-    }
-    else if (type=="normal") {
-      VaR.fm[i] <- beta.star[i,] %*% MU + 
-        sqrt(beta.star[i,,drop=F] %*% factor.star.cov %*% t(beta.star[i,,drop=F]))*qnorm(1-p)
-    }
-    # index of VaR exceedances
-    idx.exceed[[i]] <- which(R.xts <= VaR.fm[i])
-    # number of VaR exceedances
-    n.exceed[i] <- length(idx.exceed[[i]])
     
     if (type=="np") {
+      # get VaR for asset i
+      VaR.fm[i] <- quantile(R.xts, probs=1-p, na.rm=TRUE, ...)
+      
       # get F.star data object
       factor.star <- merge(factors.xts, resid.xts[,i])
       colnames(factor.star)[dim(factor.star)[2]] <- "residual"
@@ -194,11 +178,21 @@
       k.weight <- as.vector(1 - abs(R.xts - VaR.fm[i]) / eps)
       k.weight[k.weight<0] <- 0
       mVaR[i,] <- colMeans(factor.star*k.weight, na.rm =TRUE)
-    } 
-    else if (type=="normal")  {
+      
+    } else if (type=="normal")  {
+      # get VaR for asset i
+      VaR.fm[i] <- beta.star[i,] %*% MU + 
+        sqrt(beta.star[i,,drop=F] %*% factor.star.cov %*% t(beta.star[i,,drop=F]))*qnorm(1-p)
+      
+      # compute marginal VaR
       mVaR[i,] <- t(MU) + SIGB[i,] * qnorm(1-p)/sd(R.xts, na.rm=TRUE)
     }
     
+    # index of VaR exceedances
+    idx.exceed[[i]] <- which(R.xts <= VaR.fm[i])
+    # number of VaR exceedances
+    n.exceed[i] <- length(idx.exceed[[i]])
+    
     # correction factor to ensure that sum(cVaR) = asset VaR
     cf <- as.numeric( VaR.fm[i] / sum(mVaR[i,]*beta.star[i,], na.rm=TRUE) )
     
@@ -221,7 +215,6 @@
 
 fmVaRDecomp.sfm <- function(object, factor.cov, p=0.95, type=c("np","normal"), 
                             use="pairwise.complete.obs", ...) {
-  
   # set default for type
   type = type[1]
   if (!(type %in% c("np","normal"))) {
@@ -249,17 +242,14 @@
              compatible with the number of factors in the fitSfm object")
       }
     }
-    
     # get cov(F.star): (K+1) x (K+1)
     K <- object$k
     factor.star.cov <- diag(K+1)
     factor.star.cov[1:K, 1:K] <- factor.cov
     colnames(factor.star.cov) <- c(colnames(factor.cov),"residuals")
     rownames(factor.star.cov) <- c(colnames(factor.cov),"residuals")
-    
     # factor expected returns
     MU <- c(colMeans(factors.xts, na.rm=TRUE), 0)
-    
     # SIGMA*Beta to compute normal mVaR
     SIGB <- beta.star %*% factor.star.cov
   }
@@ -280,20 +270,11 @@
   for (i in object$asset.names) {
     # return data for asset i
     R.xts <- object$data[,i]
-    # get VaR for asset i
-    if (type=="np") {
-      VaR.fm[i] <- quantile(R.xts, probs=1-p, na.rm=TRUE)
-    }
-    else if (type=="normal") {
-      VaR.fm[i] <- beta.star[i,] %*% MU + 
-        sqrt(beta.star[i,,drop=F] %*% factor.star.cov %*% t(beta.star[i,,drop=F]))*qnorm(1-p)
-    }
-    # index of VaR exceedances
-    idx.exceed[[i]] <- which(R.xts <= VaR.fm[i])
-    # number of VaR exceedances
-    n.exceed[i] <- length(idx.exceed[[i]])
     
     if (type=="np") {
+      # get VaR for asset i
+      VaR.fm[i] <- quantile(R.xts, probs=1-p, na.rm=TRUE, ...)
+      
       # get F.star data object
       factor.star <- merge(factors.xts, resid.xts[,i])
       colnames(factor.star)[dim(factor.star)[2]] <- "residual"
@@ -308,11 +289,21 @@
       k.weight <- as.vector(1 - abs(R.xts - VaR.fm[i]) / eps)
       k.weight[k.weight<0] <- 0
       mVaR[i,] <- colMeans(factor.star*k.weight, na.rm =TRUE)
-    } 
-    else if (type=="normal")  {
+      
+    } else if (type=="normal")  {
+      # get VaR for asset i
+      VaR.fm[i] <- beta.star[i,] %*% MU + 
+        sqrt(beta.star[i,,drop=F] %*% factor.star.cov %*% t(beta.star[i,,drop=F]))*qnorm(1-p)
+      
+      # compute marginal VaR
       mVaR[i,] <- t(MU) + SIGB[i,] * qnorm(1-p)/sd(R.xts, na.rm=TRUE)
     }
     
+    # index of VaR exceedances
+    idx.exceed[[i]] <- which(R.xts <= VaR.fm[i])
+    # number of VaR exceedances
+    n.exceed[i] <- length(idx.exceed[[i]])
+    
     # correction factor to ensure that sum(cVaR) = asset VaR
     cf <- as.numeric( VaR.fm[i] / sum(mVaR[i,]*beta.star[i,], na.rm=TRUE) )
     
@@ -335,7 +326,6 @@
 
 fmVaRDecomp.ffm <- function(object, factor.cov, p=0.95, type=c("np","normal"), 
                             use="pairwise.complete.obs", ...) {
-  
   # set default for type
   type = type[1]
   if (!(type %in% c("np","normal"))) {
@@ -363,17 +353,14 @@
              compatible with the number of factors in the fitSfm object")
       }
     }
-    
     # get cov(F.star): (K+1) x (K+1)
     K <- ncol(object$beta)
     factor.star.cov <- diag(K+1)
     factor.star.cov[1:K, 1:K] <- factor.cov
     colnames(factor.star.cov) <- c(colnames(factor.cov),"residuals")
     rownames(factor.star.cov) <- c(colnames(factor.cov),"residuals")
-    
     # factor expected returns
     MU <- c(colMeans(factors.xts, na.rm=TRUE), 0)
-    
     # SIGMA*Beta to compute normal mVaR
     SIGB <-  beta.star %*% factor.star.cov
   }
@@ -396,20 +383,11 @@
     subrows <- which(object$data[[object$asset.var]]==i)
     R.xts <- as.xts(object$data[subrows,object$ret.var], 
                     as.Date(object$data[subrows,object$date.var]))
-    # get VaR for asset i
-    if (type=="np") {
-      VaR.fm[i] <- quantile(R.xts, probs=1-p, na.rm=TRUE, ...)
-    }
-    else if (type=="normal") {
-      VaR.fm[i] <- beta.star[i,] %*% MU + 
-        sqrt(beta.star[i,,drop=F] %*% factor.star.cov %*% t(beta.star[i,,drop=F]))*qnorm(1-p)
-    }
-    # index of VaR exceedances
-    idx.exceed[[i]] <- which(R.xts <= VaR.fm[i])
-    # number of VaR exceedances
-    n.exceed[i] <- length(idx.exceed[[i]])
     
     if (type=="np") {
+      # get VaR for asset i
+      VaR.fm[i] <- quantile(R.xts, probs=1-p, na.rm=TRUE, ...)
+      
       # get F.star data object
       factor.star <- merge(factors.xts, resid.xts[,i])
       colnames(factor.star)[dim(factor.star)[2]] <- "residual"
@@ -424,11 +402,21 @@
       k.weight <- as.vector(1 - abs(R.xts - VaR.fm[i]) / eps)
       k.weight[k.weight<0] <- 0
       mVaR[i,] <- colMeans(factor.star*k.weight, na.rm =TRUE)
-    } 
-    else if (type=="normal")  {
+      
+    } else if (type=="normal")  {
+      # get VaR for asset i
+      VaR.fm[i] <- beta.star[i,] %*% MU + 
+        sqrt(beta.star[i,,drop=F] %*% factor.star.cov %*% t(beta.star[i,,drop=F]))*qnorm(1-p)
+      
+      # compute marginal VaR
       mVaR[i,] <- t(MU) + SIGB[i,] * qnorm(1-p)/sd(R.xts, na.rm=TRUE)
     }
     
+    # index of VaR exceedances
+    idx.exceed[[i]] <- which(R.xts <= VaR.fm[i])
+    # number of VaR exceedances
+    n.exceed[i] <- length(idx.exceed[[i]])
+    
     # correction factor to ensure that sum(cVaR) = asset VaR
     cf <- as.numeric( VaR.fm[i] / sum(mVaR[i,]*beta.star[i,], na.rm=TRUE) )
     

Modified: pkg/FactorAnalytics/man/fmEsDecomp.Rd
===================================================================
--- pkg/FactorAnalytics/man/fmEsDecomp.Rd	2016-08-05 15:21:42 UTC (rev 4026)
+++ pkg/FactorAnalytics/man/fmEsDecomp.Rd	2016-08-06 18:03:52 UTC (rev 4027)
@@ -2,32 +2,42 @@
 % Please edit documentation in R/fmEsDecomp.R
 \name{fmEsDecomp}
 \alias{fmEsDecomp}
+\alias{fmEsDecomp.ffm}
 \alias{fmEsDecomp.sfm}
 \alias{fmEsDecomp.tsfm}
 \title{Decompose ES into individual factor contributions}
 \usage{
 fmEsDecomp(object, ...)
 
-\method{fmEsDecomp}{tsfm}(object, p = 0.95, type = c("np", "normal"), ...)
+\method{fmEsDecomp}{tsfm}(object, factor.cov, p = 0.95, type = c("np",
+  "normal"), use = "pairwise.complete.obs", ...)
 
-\method{fmEsDecomp}{sfm}(object, p = 0.95, type = c("np", "normal"), ...)
+\method{fmEsDecomp}{sfm}(object, factor.cov, p = 0.95, type = c("np",
+  "normal"), use = "pairwise.complete.obs", ...)
+
+\method{fmEsDecomp}{ffm}(object, factor.cov, p = 0.95, type = c("np",
+  "normal"), use = "pairwise.complete.obs", ...)
 }
 \arguments{
 \item{object}{fit object of class \code{tsfm}, \code{sfm} or \code{ffm}.}
 
 \item{...}{other optional arguments passed to \code{\link[stats]{quantile}}.}
 
+\item{factor.cov}{optional user specified factor covariance matrix with 
+named columns; defaults to the sample covariance matrix.}
+
 \item{p}{confidence level for calculation. Default is 0.95.}
 
 \item{type}{one of "np" (non-parametric) or "normal" for calculating VaR. 
 Default is "np".}
+
+\item{use}{method for computing covariances in the presence of missing 
+values; one of "everything", "all.obs", "complete.obs", "na.or.complete", or 
+"pairwise.complete.obs". Default is "pairwise.complete.obs".}
 }
 \value{
 A list containing 
 \item{ES.fm}{length-N vector of factor model ES of N-asset returns.}
-\item{n.exceed}{length-N vector of number of observations beyond VaR for 
-each asset.}
-\item{idx.exceed}{list of numeric vector of index values of exceedances.}
 \item{mES}{N x (K+1) matrix of marginal contributions to VaR.}
 \item{cES}{N x (K+1) matrix of component contributions to VaR.}
 \item{pcES}{N x (K+1) matrix of percentage component contributions to VaR.}
@@ -52,10 +62,13 @@
 contributions to \code{ES} respectively. The marginal contribution to ES is
 defined as the expected value of \code{F.star}, conditional on the loss 
 being less than or equal to \code{VaR.fm}. This is estimated as a sample 
-average of the observations in that data window.
+average of the observations in that data window. 
+
+Refer to Eric Zivot's slides (referenced) for formulas pertaining to the 
+calculation of Normal ES (adapted from a portfolio context to factor models).
 }
 \examples{
-# Time Series Factor Model
+#' # Time Series Factor Model
 data(managers)
 fit.macro <- fitTsfm(asset.names=colnames(managers[,(1:6)]),
                      factor.names=colnames(managers[,(7:8)]), data=managers)
@@ -66,14 +79,25 @@
 # Statistical Factor Model
 data(StockReturns)
 sfm.pca.fit <- fitSfm(r.M, k=2)
-ES.decomp <- fmEsDecomp(sfm.pca.fit)
+ES.decomp <- fmEsDecomp(sfm.pca.fit, type="normal")
 ES.decomp$cES
 
+# Fundamental Factor Model
+data(Stock.df)
+exposure.vars <- c("BOOK2MARKET", "LOG.MARKETCAP")
+fit <- fitFfm(data=stock, asset.var="TICKER", ret.var="RETURN", 
+              date.var="DATE", exposure.vars=exposure.vars)
+ES.decomp <- fmEsDecomp(fit, type="normal")
+head(ES.decomp$cES)
+
 }
 \author{
 Eric Zviot, Sangeetha Srinivasan and Yi-An Chen
 }
 \references{
+Eric Zivot's slides from CFRM 546: Estimating risk measures: Portfolio of 
+Assets, April 28, 2015.
+
 Epperlein, E., & Smillie, A. (2006). Portfolio risk analysis Cracking VAR 
 with kernels. RISK-LONDON-RISK MAGAZINE LIMITED-, 19(8), 70.
 

Modified: pkg/FactorAnalytics/man/fmVaRDecomp.Rd
===================================================================
--- pkg/FactorAnalytics/man/fmVaRDecomp.Rd	2016-08-05 15:21:42 UTC (rev 4026)
+++ pkg/FactorAnalytics/man/fmVaRDecomp.Rd	2016-08-06 18:03:52 UTC (rev 4027)
@@ -74,7 +74,6 @@
 data(managers)
 fit.macro <- fitTsfm(asset.names=colnames(managers[,(1:6)]),
                      factor.names=colnames(managers[,(7:8)]), data=managers)
-                     
 VaR.decomp <- fmVaRDecomp(fit.macro)
 # get the component contributions
 VaR.decomp$cVaR
@@ -82,7 +81,6 @@
 # Statistical Factor Model
 data(StockReturns)
 sfm.pca.fit <- fitSfm(r.M, k=2)
-
 VaR.decomp <- fmVaRDecomp(sfm.pca.fit, type="normal")
 VaR.decomp$cVaR
 
@@ -91,7 +89,6 @@
 exposure.vars <- c("BOOK2MARKET", "LOG.MARKETCAP")
 fit <- fitFfm(data=stock, asset.var="TICKER", ret.var="RETURN", 
               date.var="DATE", exposure.vars=exposure.vars)
-
 VaR.decomp <- fmVaRDecomp(fit, type="normal")
 VaR.decomp$cVaR
 



More information about the Returnanalytics-commits mailing list